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I.  INTRODUCTION 


This  final  report  summarizes  the  progress  made  on  research  related  to  the  problem  of 
"Pulse  Propagation  in  High  Speed  Digital  Circuits"  covering  the  period  of  September  1,  1985  to 
March  31,  1989. 


II.  PERSONNEL 

Dr.  R.  Mittra,  Professor  of  Electrical  Engineering 

Dr.  C.  Ctian,  Visiring  Assistant  Professor 

Mr.  J.  Schutt-Aine,  Visiting  Assistant  Professor 

Dr.  Z.  Pantic,  Research  Associate 

Mr.  A.  Ali,  Air  Force  Fellow 

Mr.  Paul  Aoyagi,  Graduate  Teaching  Assistant 

.Mr.  T.  Blazeck,  Graduate  Research  Assistant 

Mr.  S.  Castillo,  Graduate  Research  Assistant 

Mr.  Woo-Sung  Chi,  Graduate  Teaching  Assistant 

Mr.  U.  Feldman,  Graduate  Research  Assistant 

Mr.  P.  Harms,  Graduate  Research  Assistant 

Mr.  A.  Khebir,  Graduate  Research  Assistant 

Mr.  A.  Kouki,  Graduate  Student 

Mr.  M.  Mehalic,  Air  Force  Fellow 

-Vlr.  G.  Salo,  Graduate  Research  Assistant 

Mr.  J.  Sutton,  Graduate  Research  Assistant 

."Vlr.  G.  Wilkins,  Graduate  Student  Fellow 


HI.  TECHNICAL 

A.  QUASISTATIC  AND  FREQUENCY-DEPENDENT  ANALYSIS  OF  DIFFERENT 
TRANSMISSION-LINE  CONnGURATIONS 


We  have  developed  techniques  for  the  computation  of  the  characteristic  impedance  as  well 
as  capacitance  and  inductance  matrices  for  a  number  of  different  types  of  transmission-line 
constructions,  e.g.,  stripline,  microstrip  line,  buried  microstrip  line  and  orthogonal  line.  The 
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numerical  techniques  we  have  employed  employ  the  FFT  algorithm  combined  with  the  iterative 
techniques  and  are  considerably  more  efficient  than  those  based  upon  the  conventional  moment 
method  that  require  the  inversion  of  large  matrices.  The  algorithms  have  been  extended  to  the 
frequency-dependent  case  by  using  full-field  solution  of  Maxwell's  equations,  rather  than  the 
quasi- static  approximation. 

We  have  also  investigated  two  experimental  techniques  for  measuring  the  high  frequency 
characteristics  of  planar  transmission  lines  which  may,  in  general,  be  lossy.  In  this  approach,  the 
line  parameters,  viz.,  R,  L,  G  and  C,are  obtained  by  using  an  algorithm  for  reducing  these 
parameters  from  the  S-parameter  measurements  carried  out  with  a  Network  Analyzer. 

B.  HNITE  ELEMENT  ANALYSIS  OF  TRANSMISSION  LINES 

We  have  developed  a  general  purpose,  finite  element  (FEM)  analysis  which  is  capable  of 
handling  transmission  lines  of  arbitrary  cross-section  and  with  arbitrary  fillings.  The  FEM  program 
requires  a  mesh  generation  algorithm.  A  simple  mesh  program,  that  is  portable  in  nature,  was 
developed  by  us  for  use  with  the  program.  However,  we  have  found  that  it  is  more  expedient  to 
acquire  a  commercial  version  of  the  mesh  generation  program,  e.g.,  PATRAN,  and  combine  it 
with  our  FEM  programs  for  solving  either  Laplace's  or  Maxwell's  equations. 

The  FEM  approach  calls  for  a  method  for  mesh  truncation  whenever  open-region  problems 
are  analyzed  using  this  method.  We  have  developed  an  absorbing  boundary  condition  for  mesh 
truncation  and  have  applied  it  to  single  and  multi-conductor  transmission  lines.  The  approach  is 
described  in  Attachment  A. 

C.  PULSE  PROPAGATION  AND  CROSSTALK  IN  MULTICONDUCTOR  LINES 

WITH  LOSSES 

We  have  generalized  the  analysis  of  lossless  n-lines,  terminated  with  digital  devices,  to  the 
lossy  case  using  two  different  approaches.  In  the  first  approach,  the  frequency-domain  Gicen'.s 
function  method  is  used  in  the  transmission  line  part  of  the  circuit,  which  is  linear,  in  conjunction 
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with  a  nonlinear  analysis  of  the  terminating  circuits,  that  include  digital  devices,  directly  in  the  time 
domain.  The  line  is  characterized  by  per  unit  values  of  resistance  (r),  inductance  (1),  conductance 
(g)  and  capacitance  (c).  If  desired,  r  can  be  made  to  vary  with  frequency  to  model  the  skin  effect 
and  g  can  also  be  a  variable  of  frequency.  The  terminal  networks  that  model  the  nonlinear  digital 
devices  consist  of  capacitors  in  parallel  with  nonlinear  voltage-controlled  resistors. 

We  have  prepared  a  report  (see  Attachment  B)  describing  this  analysis  and  have  appended  it 
herewith. 

The  second  approach  we  have  developed  for  the  time-domain  simulation  of  n-line 
structures  terminated  by  logic  gate  inverters,  is  based  on  a  state- space  approach  employed  in 
conjunction  with  the  Scattering  Matrix  representation  of  the  coupled  n-lines.  It  is  designed  to  treat 
both  the  lossless  and  lossy  line  cases  and  has  been  described  in  a  recent  publication  a  copy  of 
which  is  appended  in  Attachment  C. 

D.  CONNECTOR,  INTERCONNECT  AND  DISCONTINUITY  MODELING 

The  finite  element  method,  as  well  as  other  analytical  and  numerical  techniques,  are  being 
investigated  with  the  objective  of  modeling  connectors,  interconnects  and  discontinuities 
introduced  by  a  change  in  the  conductor  width  or  the  presence  of  a  stub,  etc.  In  general,  the 
connector,  the  via,  and  similar  other  discontinuity  problems  must  be  modeled  by  using  a  three- 
dimensional  mesh  generation  program;  hence,  such  modeling  is  computer  intensive.  The 
connector  problem  has  been  investigated  using  the  full  3-D  mesh  and  with  an  approximate  2-D 
model.  It  has  been  shown  that  in  some  cases  it  is  possible  to  use  the  2-D  approximation,  though 
this  is  not  universally  true.  A  comprehensive  report  on  the  connector  that  describes  the  analysis, 
numerical  results  and  the  experimental  verification  of  these  results  has  been  prepared  and  a  paper 
summarizing  the  research  has  been  submitted  for  publication. 

An  alternative  approach,  based  on  the  method  of  moments  (MoM)  has  been  employed  for 
both  the  bend  and  the  via  problems.  Numerical  results  for  the  inductance  and  capacitance  matrices 
have  been  obtained  for  the  single  and  multiple  line  bends  in  microstrip  lines.  The  modeling  of  of 
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the  via  problem  is  being  continued. 

Another  discontinuity  problem  being  investigated  in  this  effort  is  the  tapered-line  problem, 
which  arises  in  the  modeling  of  interconnects,  and  has  been  analyzed  using  a  perturbational 
approach.  Both  single  and  multiple  lines  with  arbitrary  tapers  have  been  investigated  using  the 
scattering  matrix  approach.  Experimental  verification  of  the  theoretical  results  is  currently  being 
carried  out. 
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APPENDIX  A 


AN  ABSORBING  BOUNDARY  CONDITION 
FOR  QUASI-TEM  ANALYSIS  OF  MICROWAVE 
TRANSMISSION  LINES  VIA  THE  FINITE  ELEMENT 

METHOD 


A.  Khebir,  A.  B.  Kouki,  and  R.  Mittra 
Electromagnetic  Communication  Laboratory 
Electrical  &  Computer  Engineering  Department 
University  of  Illinois 
Urbana,  IL  61801 


Abstract 

Open  microwave  transmission  line  structures  are  analyzed  in  this  paper  in  the  quasi- 
TEM  regime  using  the  finite  element  method.  An  absorbing  boundary  condition(ABC)  for  an 
arbitrary  outer  boundary  is  introduced  for  the  purpose  of  truncating  the  unbounded  region 
surrounding  the  transmission  line  in  an  efficient  manner.  The  application  of  the  boundary 
condition  is  illustrated  for  three  different  microstrip  line  configurations,  viz.,  a  single  line,  two 
coupled  lines,  and  a  six-conductor  line.  Numerical  results  are  compared  with  those  published 
elsewhere  and  good  agreement  is  found. 
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I.  Introduction 

Microwave  transmission  lines  have  been  investigated  by  many  researchers  who  have 
employed  a  variety  of  methods  to  study  the  problem  of  computing  the  characteristic  impedance 
and  propagation  constant  along  these  lines.  Some  of  these  techniques  include  the  Fourier 
transform  method  [1-2],  variational  method  [3-4],  spectral  domain  method  [5-7],  Green's 
function  technique  [8-13],  conformal  mapping  [14-16],  boundary  element  method  [17-18],  and 
finite-element  method  [19-20].  All  but  the  last  two  approaches  mentioned  above  are  limited  to 
thin  strips  and/or  to  structures  containing  dielectrics  with  planar  interfaces.  Although  the  finite 
element  method  method  (FEM)  is  very  general,  and  can  handle  any  arbitrary  configuration  of 
conductors  and  dielectrics,  it  must  deal  with  the  practical  problems  of  mesh  truncation  and  the 
need  for  a  large  number  of  mesh  nodes  when  applied  to  an  open  region  problem.  One  approach 
to  circumventing  this  difficulty  is  to  truncate  the  mesh  by  introducing  a  fictitious  conducting 
enclosure[l9,  21],  This  approach  yields  satisfactory  results  only  if  the  actual  field  decays 
sufficiently  well  as  it  reaches  the  outer  boundary.  Typically,  this  requires  one  to  recede  the 
outer  boundary  far  away  from  the  structure  in  order  to  achieve  acceptable  accuracy  and,  this,  in 
turn,  results  in  a  large  mesh.  An  alternative  approach  is  to  use  "infinite"  elements  [20],  that 
extend  to  infinity,  and  cover  the  region  outside  of  a  fictitious  boundary  surrounding  the 
structure.  Although  superior  to  the  artificial  p.e.c  boundary  method,  this  approach  nonetheless 
has  its  own  drawbacks.  First,  the  infinite  elements  require  special  care  during  the  filling  of  the 
FEM  matrix.  Second,  one  needs  to  assume  a  certain  asymptotic  behavior  of  the  field  within  the 
infinite  elements,  and  this  behavior  may  not  be  convenient  to  obtain. 

In  this  paper,  we  introduce  an  absorbing  boundary  condition(ABC)  which  provides  us 
with  an  efficient  means  for  dealing  with  the  open  region  problems  in  the  quasi-static  regime. 
This  absorbing  boundary  condition  does  not  suffer  from  the  complications  associated  with  the 
infinite  elements,  and  yet  enables  us  to  bring  the  outer  boundary  much  closer  to  the  structure 
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than  would  be  possible  with  the  p.e.c.  artificial  boundary.  Furthermore,  unlike  many  of  the 
available  ABC's  that  are  restricted  to  separable  outer  boundaries,  the  one  presented  in  this 
paper  is  useful  for  an  arbitrarily-shaped  outer  boundary.  We  will  demonstrate  the  versatility  of 
this  new  absorbing  boundary  condition  by  considering  the  examples  of  one,  two,  and  six 
conductor  microstrip  lines. 

II.  Derivation  of  the  Absorbing  Boundary  Condition 

Figure  1  depicts  the  geometry  of  an  open  region  problem  consisting  of  N  arbitrarily- 
shaped  conductors  embedded  in  a  multilayered  medium  above  a  ground  plane.  Let  fl-r  denote 

the  region  exterior  to  the  conductors.  For  the  finite  mathematics  techniques  [22-23],  the 
unbounded  outer  region  Qj  must  be  truncated  and  enclosed  with  an  outer  boundary  I2.  The 

strategy  for  deriving  an  absorbing  boundary  condition  is  to  find  an  operator  which,  when 
applied  on  the  outer  boundary,  would  mimic  the  asymptotic  behavior  of  the  wave  function  at 
,V;ity  and  would  thus  yield  accurate  results  in  the  interior  region  without  the  need  of  an 
exorbitantly  large  number  of  mesh  points.  The  problem  at  hand  can  be  described  in  terms  of 
the  following  set  of  equations: 

V^u  =  0  in  (1) 

u  =  gj  on  the  i'^  conductor  (2) 


B„u  =  0  on  r, 
m  z 


where  u  is  the  potential  function  and  Bm  is  the  m'*’  order  absorbing  boundary  condition 
operator. 


An  asymptotic  form  of  the  solution,  valid  for  large  p  can  be  written  as: 


a^ 

u  (p,(l))  =2^  —  cos  n<l) 

n=l  P" 


(4) 


Equation  (4)  can  be  used  to  obtain  an  absorbing  boundary  condition  on  an  artificial  boundary. 
Manipulating  (4)  one  can  obtain  the  first-order  absorbing  boundary  condition  operator,  which 
reads 


Similarly,  the  second-order  absorbing  boundary  condition  operator  takes  the  form 

As  will  be  seen  later,  the  finite  element  formulation  of  the  problem  involves  an  integral  over  the 
outer  boundary,  I2,  and  the  normal  derivative  of  u  appears  in  this  integrand.  Hence,  for  our 

purposes,  it  it  is  more  desirable  to  find  an  asymptotic  representation  for  the  normal  derivative 
of  u  rather  than  make  direct  use  of  the  operator  B2.  For  a  circular  outer  boundary,  the  normal 
derivative  will  simply  be  the  radial  one  which  can  readily  be  obtained  from  (6).  Mittra  and 
Ramahi  [30]  have  presented  a  procedure  for  deriving  an  absorbing  boundary  condition  operator 
for  the  dynamic  case  in  terms  of  the  field  and  its  second-order  angular  derivative  that  can  be 
conveniently  adapted  to  the  quasi-static  case  of  interest  here.  We  begin  by  expressing  the  radial 
derivative  of  u  in  a  manner  similar  to  the  one  described  in  [30].  Next,  we  use  (1)  to  exchange 
the  second-order  derivative  in  p,  Upp,  with  the  second-order  angular  derivative,  U(jK)),  to  obtain 
the  desired  absorbing  boundary  condition: 

Up  =  a(p)  u  -f-  p(p)  u^^  (7) 


where  a(p)  and  P(p)  are  given  by 
a(p)  =  ~ 

P(P)  =  ^ 


(8) 


(9) 
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For  a  general,  noncircular  outer  boundary  it  is  necessary  to  generalize  (7)  and  derive  an 
expression  for  the  normal  derivative  operator  Un  in  the  local  coordinate  system  (t,n)  where  t  and 
n  are  tangent  and  normal  to  the  boundaiy,  respectively.  An  approximate  expression  for  this 
derivative  can  be  obtained  by  following  the  procedure  described  in  [24]  and  is  given  by 

|^  =  au+ Y  Uj+p  u„  (10) 


where 


—  f  2  \ 

a  =  (xq  sinOg  -  Yq  cosOq) - - 

V  3p  > 


1  2  2 

_  -t  (xq  sinOg  -yg  cosOg)  +-jsin20g  (yg  -Xg)  +Xgyo  cos29g 
7= - ; - 


—  3/^1' 

P  =  (Xg  sinGg  -yg  cosOg)  - r 

Up  . 


+x^  +  yg  +  2t  (xgCosOg  +ygSin8o) 


(11) 

(12) 

(13) 

(14) 


and  00,  XQ,  yo.  and  t  are  as  shown  in  Fig.  2. 


III.  Finite  Element  Implementation  of  the  Absorbing  Boundary  Condition 


The  problem  at  hand  is  to  solve  for  the  potential  u  satisfying  the  Laplace's  equation; 


V-(eVu)  =  0 


(15) 


Multiplying  (15)  by  a  testing  function  v  and  integrating  over  the  domain  of  the  problem 
Qj,  we  obtain 


[  V  V-(eVu)  ds  =  0 


(16) 
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From  Green's  identity  we  have 


f  V  V  (eVu)  ds  =  -  f  eVu-Vv  ds  +  f  v  dt 
J^r  Jfi,  /  dn 


Inserting  the  above  in  (16),  we  get 

9u 


I  eVu-Vv  ds  = 

Jo, 


r  ou 


(17) 


(18) 


In  the  finite  element  formulation,  one  sets  up  a  mesh  in  the  region  typically  using 
triangular  elements.  The  edges  of  the  outermost  elements  prescribe  r2.  Hence,  considering 
one  element  at  a  time,  the  absorbing  boundary  condition  given  in  (10)  may  be  incorporated  into 
(18)  to  yield 


eVu-Vv  ds  =  I  V  e  ( a  u  +  Y  Uj  +  (5  u^)  dt 


(19) 


Since  e  is  constant  over  each  element,  we  can  inte^te  (19)  by  parts  to  obtain 


J  eVu'Vv  ds  =  J  e  (a  vu  +  Y  vUj  -  p  VjUj  -  vUj)  dt 


(20) 


where 

^  =  (xoSin0o-yQCOs0o)^(t  +  XoCos0o+yoSin0Q)r — (21) 

V  3p  ; 

The  form  given  in  (20)  is  well-suited  for  numerical  implementation  for  any  arbitrary 
conductor  configuration  enclosed  by  any  arbitrarily-shaped  outer  boundary. 
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IV.  Numerical  Results 

a.  One  Conductor 

The  microstrip  line,  shown  in  Fig.  3,  is  enclosed  by  a  rectangular  outer  boundary.  We 
first  solved  the  potential  problem  by  applying  the  absorbing  boundary  condition,  given  in 
equation  (7),  on  a  rectangular  outer  boundary.  Next,  we  introduced  a  perfectly  electric 
conducting  shield  at  the  outer  boundary  and  solved  the  problem  once  again  using  the  same 
mesh.  As  Table  1  indicates,  the  relative  error  between  the  absorbing  boundary  condition  and 
the  published  results  [12,  26,  27]  is  between  0.053  and  2.56  percent,  whereas  the  error 
between  the  shield  and  the  published  results  is  between  14.55  and  30.04  percent. 

b.  Two  conductors: 

Two  coupled  microstrips,  shown  in  Fig.  4,  are  enclosed  by  a  rectangular  outer 
boundary.  Table  2  presents  some  results  for  the  same  problem  that  have  been  published 
elsewhere(12,  28, 29],  together  with  those  obtained  by  using  a  p.e.c.  shield  and  the  absorbing 
boundary  condition  in  (10).  The  relarive  error  in  the  self-terms  is  0.27  percent  for  the 
absorbing  boundary  condition  and  18.31  percent  for  the  shield.  Likewise,  the  error  in  the 
mutual  terms  is  5.2  percent  for  the  absorbing  boundary  condition  and  44.59  percent  for  the 
shield. 

c. Six  conductors: 

The  six  conductor  system,  shown  in  Fig.  5,  is  enclosed,  once  again,  by  a  rectangular 
outer  boundary.  To  the  best  of  our  knowledge,  there  are  no  published  results  for  this 
configuration.  However,  we  have  compared  our  results  with  those  derived  by  using  tne 
computer  program  developed  by  Chan  [25],  which  uses  an  integral  equation  formulation  and 
an  iterative  method  of  solution.  As  Table  3  indicates,  the  relative  error  for  the  capacitance 
matrix  is  between  0.84  and  14.52  percent  for  the  absorbing  boundary  condition  and  between 
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3.96  and  74.10  percent  for  the  p.e.c  shield  when  compared  to  the  results  obtained  with  the 
iterative  solution.  For  the  inductance  matrix,  the  error  is  between  0.02  and  11.46  percent  for 
the  absorbing  boundary  condition  and  between  39.74  and  96.72  percent  for  the  p.e.c  shield. 
In  addition,  the  relative  error  between  the  absorbing  boundary  condition  and  the  shield  reached 
69.69  percent  in  the  capacitance  matrix  and  97.01  percent  in  the  inductance  matrix. 

In  all  of  the  three  numerical  examples  considered  thus  far,  we  have  chosen  a  rectangular 
outer  boundary  because  it  facilitates  the  meshing  procedure  and  because  it  is  conformal  to  the 
structures  considered.  This  meshing  is  done  in  a  manner  such  that  none  of  the  triangular 
elements  have  more  than  one  edge  on  the  outer  boundary.  Since  the  finite  element  scheme 
considers  one  element  at  a  time,  the  problem  of  the  undefined  normal  at  the  rectangular  comers 
is  circumvented  when  the  procedure  described  above  is  followed. 

In  order  to  illustrate  the  fact  that  the  formulation  described  in  this  paper  is  also 
applicable  to  an  arbitrary  outer  boundary,  we  reconsider  the  two-conductor  example  described 
earlier.  Figure  6  illustrates  the  coupled  microstrips  problem  with  an  arbitrary  outer  boundary. 
Table  4  shows  that  the  results  obtained  for  this  case  are  almost  identical  to  those  derived  with 
the  rectangular  boundary. 

Finally,  before  closing,  we  point  out  that  no  special  treatment  is  needed  at  the  dielectric 
interfaces  because  in  the  finite  element  formulation  the  medium  is  modeled  as  being 
homogeneous  within  each  element  and,  consequently,  the  line  integral  in  equation  (20)  is 
always  confined  to  within  a  homogeneous  region  inside  the  element. 

V.  Conclusions 

We  have  shown  how  an  absorbing  boundary  condition  for  quasi-static  fields  can  be 
applied  to  a  potential  field  at  distances  quite  close  to  a  transmission  line  configuration  to  derive 


an  FEM-based  solution  in  a  numerically  efficient  manner.  It  is  evident  from  the  numerical 
results  that  the  absorbing  boundary  condition  consistently  yields  more  accurate  results  than 
those  obtainable  with  a  perfectly  conducting  shield  placed  at  the  same  location.  However,  in 
some  situations,  the  accuracy  obtained  with  the  approximate  ABC  presented  in  this  paper  may 
still  not  be  adequate,  as  for  instance  in  the  off-diagonal  terms  of  the  capacitance  matrix  in  the 
six-conductor  configuration.  This  aspect  of  the  problem  is  being  investigated  further  by 
attempting  to  derive  Improved  versions  of  the  absorbing  boundary  condition.  Also  planned  is 
an  extension  to  the  full-wave  FEM  formulation  where  a  dynamic  version  of  the  absorbing 
boundary  condition  is  needed. 
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Tabic  3a:  Capacitance  matrix  for  the  six-condutor  structure  of  Hg.  5. 


Iterative 

[25] 

ABC 

Shield 

Error  Error 

ABC-[25]  Shield-[25] 

C(l,l) 

0.668  X  10-10 

0.686  X  10-10 

0.848  X  10-10 

2.620% 

26.83% 

C(1.2) 

-0.279  X  10-10 

-0.315  X  10-10 

-0.264  X  10-10 

13.05% 

5.340% 

C(1.3) 

-0.549  X  10-11 

-0.600  X  10-11 

-0.371  X  10-11 

9.240% 

32.56% 

C(l,4) 

-0.208  X  10-11 

-0.225  X  10-11 

-0.117  X  10-11 

8.320% 

43.71% 

C(l,5) 

-0.999  X  10-12 

-0.101  X  10-11 

-0.456  X  10-12 

0.840% 

54.30% 

C(l,6) 

-0.704  X  10-12 

-0.602  X  10-12 

-0.182  X  10-12 

14.52% 

74.10% 

C(2,l) 

-0.279  X  10-10 

-0.315  X  10-10 

-0.264  X  10-10 

13.05% 

5.34% 

C(1J1) 

0.789  X  10-10 

0.848  X  10-10 

0.876  X  10-10 

7.480% 

11.03% 

C(2.3) 

-0.256  X  10-10 

-0.284  X  10-10 

-0.266  X  10-10 

11.10% 

3.960% 

C(2,4) 

-0.465  X  10-11 

-0.487  X  10-11 

-0.385  X  10-11 

4.740% 

17.23% 

C(2,5) 

-0.173  X  10-11 

-0.181  X  10-11 

-0.122  X  10-11 

4.610% 

29.03% 

C(2,6) 

-0.999  X  10-12 

-0.101  X  10-11 

-0.456  X  10-12 

0.920% 

54.30% 

C(3,l) 

-0.549  X  10-11 

-0.600  X  10-11 

-0.371  X  10-11 

9.240% 

32.56% 

C(3,2) 

-0.256  X  10-10 

-0.284  X  10-10 

-0.266  X  10-10 

11.10% 

3.960% 

C(3.3) 

0.794  X  10-10 

0.855  X  10-10 

0.874  X  10-10 

7.680% 

10.14% 

C(3.4) 

-0.254  X  10-10 

-0.282  X  10-10 

-0.267  X  10-10 

10.91% 

4.870% 

C(3,5) 

-0.465  X  10-11 

-0.487  X  10-11 

-0.385  X  10-11 

4.750% 

17.23% 

C(3.6) 

-0.208  X  10-11 

-0.226  X  10-11 

-0.117  X  10-11 

8.410% 

43.71% 

C(4,l) 

-0.208  X  10-11 

-0.225  X  10-11 

-0.117  X  10-11 

8.320% 

43.71% 

C(4,2) 

-0.465  X  10-11 

-0.487  X  10-11 

-0.385  X  10-11 

4.740% 

17.23% 

C(4.3) 

-0.254  X  10-10 

-0.282  X  10-10 

-0.267  X  10-10 

10.91% 

4.870% 

C(4,4) 

0.794  X  10-10 

0.855  X  10-10 

0.874  X  10-10 

7.680% 

10.14% 

C(4.5) 

-0.256  X  10-10 

-0.284  X  10-10 

-0.266  X  10-10 

11.10% 

3.960% 

Table  3a  (continued):  Capacitance  matrix  for  the  six-condutor  structure  of  Fig.  5. 


Iterative 

[25] 

ABC 

Shield 

Error  Error 

ABC-[25]  Shield-[25] 

C(4.6) 

-0.549  X  10-11 

-0.601  X  10-11 

-0.371  X  10-11 

9.310% 

32.56% 

C(5.1) 

-0.999x  10-12 

-0.101  X  10-11 

-0.456  X  10-12 

0.840% 

54.30% 

C(5.2) 

-0.173  X  10-11 

-0.181  X  10-11 

-0.122  X  10-11 

4.610% 

29.03% 

C(5.3) 

-0.465  X  10-11 

-0.487  X  10-11 

-0.385  X  10-11 

4.750% 

17.23% 

C(5,4) 

-0.256  X  10-10 

-0.284  X  10-10 

-0.266  X  10-10 

11.10% 

3.960% 

C(5,5) 

0.789  X  10-10 

0.848  X  10-10 

0.876  X  10-10 

7.480% 

11.03% 

C(5,6) 

-0.279  X  10-10 

-0.316  X  10-10 

-0.264  X  10-10 

13.08% 

5.340% 

C(6,l) 

-0.704  X  10-12 

-0.602  X  10-12 

-0.182  X  10-12 

14.52% 

74.10% 

C(6.2) 

-0.999  X  10-12 

-0.101  X  10-11 

-0.456  X  10-12 

0.920% 

54.30% 

C(6,3) 

-0.208  X  10-11 

-0.226  X  10-11 

-0.1 17  X  10-11 

8.410% 

43.71% 

C(6,4) 

-0.549  X  10-11 

-0.601  X  10-11 

-0.371  X  10-11 

9.310% 

32.56% 

C(6,5) 

-0.279  X  10-10 

-0.316  X  10-10 

-0.264  X  10-10 

13.08% 

5.340% 

C(6.6) 

0.668  X  10-10 

0.685  X  10-10 

0.848  X  10-10 

2.550% 

26.83% 

Table  3b:  Inductance  matrix  for  the  six-condutor  structure  of  Fig.  5. 


Iterative 

[25] 

ABC 

Shield 

Error  Error 

ABC-[25]  ShieId-[25] 

L(l,l) 

0.843  X  10-6 

0.884  X  10-6 

0.475  X  10-6 

4.830% 

43.60% 

L(1.2) 

0.439  X  10-6 

0.486  X  10-6 

0.136  X  10-6 

10.66% 

68.99% 

L(1.3) 

0.309  X  10-6 

0.344  X  10-6 

0.542  X  10-2 

11.26% 

82.44% 

L(1.4) 

0.235  X  10-6 

0.260  X  10-6 

0.242  X  10-7 

10.72% 

89.68% 

L(1.5) 

0.186  X  10-6 

0.205  X  10-6 

0.113  X  10-7 

9.830% 

93.92% 

L(1.6) 

0.153  X  10-6 

0.167  X  10-6 

0.501  X  10-8 

9.640% 

96.72% 

L(2.1) 

0.439  X  10-6 

0.486  X  10-6 

0.136  X  10-6 

10.66% 

68.99% 

L(2.2) 

0.837  X  10-6 

0.836  X  10-6 

0.497  X  10-6 

0.070% 

40.52% 

L(2.3) 

0.436  X  10-6 

0.455  X  10-6 

0.147  X  10-6 

4.300% 

66.18% 

LC2,4) 

0.307  X  10-6 

0.325  X  10-6 

0.598  X  10-7 

5.950% 

80.55% 

L(2,5) 

0.234  X  10-6 

0.252  X  10-6 

0.265  X  10-7 

7.650% 

88.69% 

L(2.6) 

0.186  X  10-6 

0.205  X  10-6 

0.113  X  10-7 

9.960% 

93.91% 

L(3.1) 

0.309  X  10-6 

0.344  X  10-6 

0.542  X  10-7 

11.26% 

82.44% 

L(3.2) 

0.436  X  10-6 

0.455  X  10-6 

0.147  X  10-6 

4.300% 

66.18% 

L(3.3) 

0.835  X  10-6 

0.816  X  10-6 

0.503  X  10-6 

2.260% 

39.74% 

L(3,4) 

0.436  X  10-6 

0.446  X  10-6 

0.150  X  10-6 

2.350% 

65.59% 

L(3.5) 

0.307  X  10-6 

0.325  X  10-6 

0.598  X  10-7 

6.020% 

80.55% 

L(3.6) 

0.235  X  10-6 

0.261  X  10-6 

0.242  X  10-7 

10.91% 

89.68% 

L(4.1) 

0.235  X  10-6 

0.260  X  10-6 

0.242  X  10-7 

10.72% 

89.68% 

L(4.2) 

0.307  X  10-6 

0.325  X  10-6 

0.598  X  10-7 

5.950% 

80.55% 

L(4,3) 

0.436  X  10-6 

0.446  X  10-6 

0.150  X  10-6 

2.350% 

65.59% 

L(4,4) 

0.835  X  10-6 

0.817  X  10-6 

0.503  X  10-6 

2.250% 

39.74% 

U4,5) 

0.436  X  10-6 

0.455  X  10-6 

0.147  X  10-6 

4.370% 

66.17% 

Table  3b  (continued):  Inductance  matrix  for  the  six-condutor  structure  of  Rg.  5. 


Iterative 

[25] 

ABC 

Shield 

Error  Error 

ABC-[25]  Shield-[25] 

L(4.6) 

0.309  X  10-6 

0.344  X  10-6 

0.542  X  10-7 

11.46% 

82.44% 

L(5,l) 

0.186  X  10-6 

0.205  X  10-6 

0.113  X  10-7 

9.830% 

93.92% 

L(5.2) 

0.234  X  10-6 

0.252  X  10-6 

0.265  X  10-7 

7.650% 

88.69% 

L(5.3) 

0.307  X  10-6 

0.325  X  10-6 

0.598  X  10-7 

6.020% 

80.55% 

L(5,4) 

0.436  X  10-6 

0.455  X  10-6 

0.147  X  10-6 

4.370% 

66.17% 

L(5,5) 

0.837  X  10-6 

0.836  X  10-6 

0.497  X  10-6 

0.020% 

40.52% 

L(5,6) 

0.439  X  10-6 

0.487  X  10-6 

0.136  X  10-6 

10.82% 

68.99% 

L(6,l) 

0.153  X  10-6 

0.167  X  10-6 

0.501  X  10-8 

9.640% 

96.72% 

L(6.2) 

0.186  X  10-6 

0.205  X  10-6 

0.113  X  10-7 

9.960% 

93.91% 

L(6.3) 

0.235  X  10-6 

0.261  X  10-6 

0.242  X  10-7 

10.91% 

89.68% 

L(6,4) 

0.309  X  10-6 

0.344  X  10-6 

0.542  X  10-7 

11.46% 

82.44% 

L(6,5) 

0.439  X  10-6 

0.487  X  10-6 

0.136  X  10-6 

10.82% 

68.99% 

L(6.6) 

0.843  X  10-6 

0.885  X  10-6 

0.475  X  10-6 

4.940% 

43.60% 

20 


ABC  with 
rectangular  outer 
boundary 

Reference 

[29] 

ABC  with 
arbitrary  outer 
boundary 

C(l,l) 

0.9249  X  10-10 

0.9224  X  10-10 

0.9284  X  10-10 

C(1.2) 

-0.8061x10-11 

-0.8504x10-11 

-0.8036x10-11 

C(2,l) 

-0.8061x10-11 

-0.8504x10-11 

-0.8036x10-11 

C(2.2) 

0.9249  X  10-10 

0.9224x10-10 

0.9284x10-10 

Table  4:  Capacitance  matrix  for  the  coupled  microstrips  of  Rgures  4  and  6. 
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CHAPTER  1 

INTRODUCTION 

Signal  propagation  along  multiconductor  transmission  lines  is  of  interest  in  the 
design  of  digital  computers  and  communication  networks.  With  the  increase  in  switching 
speeds  and  denser  packaging  of  digital  devices,  it  has  become  necessary  to  consider  trans¬ 
mission  line  effects  in  the  interconnection  of  high-speed  circuits.  Reflections  from  un¬ 
matched  terminations  or  discontinuities,  or  crosstalk  from  neighboring  lines,  may  degrade 
pulses  sufficiently  to  cause  false  logic  levels.  Skin-effect  and  dielectric  losses  cause  dis¬ 
persion  and  attenuation  that  can  increase  pulse  rise  and  fall  times,  and  interconnection 
lengths  introduce  signal  delays  that  can  limit  overall  system  speed.  These  are  among  the 
problems  confronting  designers  of  transmission  systems  that  propagate  pulses  with  short 
rise  and  fall  times  and  narrow  widths. 

High-speed  interconnects  can  be  metal  or  polysilicon  lines  on  integrated  circuit 
chips,  printed  wires  on  printed  circuit  boards,  or  microstrip  lines  between  boards.  In  addi¬ 
tion,  it  may  be  necessary  to  consider  the  effects  of  connectors,  bends,  junctions, 
crossovers,  and  other  discontinuities  in  the  transmission  line  structure.  The  logic  circuits 
that  terminate  the  lines  are  composed  of  highly  nonlinear  active  devices  such  as  transistors 
and  diodes.  The  result  can  be  a  very  complicated  network  consisting  of  lumped  linear  and 
nonlinear  elements  in  combination  with  lossy  coupled  transmission  lines.  The  goal  of  this 
research  is  to  develop  analytical  and  numerical  methods  for  predicting  the  transient  re¬ 
sponse  of  complex  networks  of  this  type.  Current  computer-aided  design  tools  for  nonlin¬ 
ear  circuits  have  the  ability  to  model  single  lossless  lines  or  use  RC  ladders  to  model  single 
resistive  lines.  They  cannot  accurately  portray  the  multiline  coupling  and  frequency- 
dependent  losses  that  arise  when  high-speed  signals  propagate  between  digital  circuits. 
Although  the  techniques  presented  herein  can  be  applied  to  the  transient  analysis  of  other 
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multiconductor  transmission  lines,  such  as  power  lines  under  short-circuit  or  power-on 
conditions,  the  principal  concern  here  is  with  the  high-speed  interconnect  lines. 

Early  researchers  of  multiple  parallel  transmission  lines  used  transform  or  operator 
techniques  to  arrive  at  analytical  solutions  for  lossless  or  lossy  lines  with  linear  termina¬ 
tions  [l]-[5].  Carson  and  Hoyt  [1]  formulated  the  multiconductor  transmission  line  equa¬ 
tions  in  the  frequency  domain.  Their  formulation  was  sufficiently  general  to  allow  for  the 
frequency  dependence  of  the  resistance  and  inductance  that  results  from  the  skin  and 
proximity  effects  [6]-[7].  Shortly  thereafter,  Bewley  [2]  included  the  concept  of  multive¬ 
locity  waves  in  a  comprehensive  analysis  of  lossless  and  lossy  multiconductor  transmission 
lines.  Later,  Pipes  [3]  and  Rice  [4]  used  a  matrix  formulation  to  simplify  the  mechanics  of 
solving  transmission  line  problems.  And  more  recently,  Wedepohl  [5]  investigated  the 
properties  of  the  characteristic  impedance  matrix,  modal  transformation  matrices,  and 
propagation  modes  associated  with  the  steady-state  solution  of  lossy  multiconductor  trans¬ 
mission  lines. 

The  early  sixties  saw  the  introduction  of  digital-computer  methods  to  aid  in  the 
solution  of  transmission  line  problems.  Transient  analysis  of  single  and  multiple  lossy 
lines  with  linear  loads  was  accomplished  by  using  the  Fourier  transform  to  switch  to  the 
time  domain  [8]-[10].  However,  as  previously  noted,  the  interconnections  of  digital  cir¬ 
cuits  are  themselves  best  modeled  as  multiconductor  transmission  lines,  and  digital  circuits 
add  to  the  problem  because  they  are  nonlinear.  Initial  concerns  dealt  with  signal  delays, 
reflections,  and  coupling  in  lossless  two-line  configurations  [1 1]-[14].  Experimental  and 
analytical  investigations  were  carried  out  that  increased  the  qualitative  understanding  of 
these  effects. 

Several  studies  of  lossless  n-conductor  transmission  lines  with  arbitrary  nonlinear 
terminations  were  conducted  in  the  late  sixties  and  early  seventies  [15]-[19].  Amemiya 
[15]  decoupled  the  propagating  signals  into  forward  and  backward  modal  waves  and  then 
employed  the  generalized  Thevenin  or  Norton  equivalent  circuits  associated  with  the 
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method  of  characteristics  to  model  the  multiconductor  line  at  its  endpoints.  Dommel  [16] 
included  this  procedure  in  an  electromagnetic  transients  program.  Chang  [17]  used  con¬ 
gruence  transformers  to  decouple  n-conductor  transmission  lines  into  n  single  lines. 
Recently,  Tripathi  and  Rettig  [20]  implemented  this  method  in  the  popular  SPICE  program. 
Ho  [18]  and  Marx  [19]  analyzed  the  properties  of  the  propagation  modes,  equivalent  cir¬ 
cuits,  and  characteristic  resistance  and  conductance  matrices  for  lossless  multiconductor 
transmission  lines.  Ho  also  presented  some  numerical  techniques  for  variable  time  step 
control  during  computer  solution  of  the  transmission  line  equations. 

Lossy  multiple  parallel  transmission  lines  with  nonlinear  terminations  first  received 
attention  in  the  late  sixties  [21]-[24].  Silverberg  and  Wing  [21]  outlined  a  methodology  for 
solving  circuits  containing  lumped  nonlinear  elements  and  distributed  linear  elements. 
Their  idea  consisted  of  first  characterizing  the  distributed  part  of  the  network  in  the  fre¬ 
quency  domain,  then  converting  this  characterization  to  the  time  domain,  and  finally  solv¬ 
ing  the  whole  system  in  a  time-stepping  manner  by  matching  boundary  conditions  at  the 
interface.  Short  circuit  admittance  parameters  (Y  parameters)  were  used  to  characterize  the 
distributed  elements.  In  the  time  domain,  these  become  impulse  responses  that  must  be 
convolved  with  the  line  voltages.  Budner  [22]  used  Y  parameters  to  show  results  for  a 
symmetric  two-line  configuration,  and  Snelson  [23]  modified  the  method  of  characteristics 
to  get  impulse  responses  of  shorter  duration  than  those  for  the  Y  parameters.  Shorter 
impulse  responses  reduce  the  problem  of  aliasing  encountered  when  taking  the  inverse 
transform  and  lessen  the  number  of  calculations  that  need  to  be  computed  in  the  time 
domain.  Meyer  and  Dommel  [24]  made  several  numerical  improvements  to  Snelson’s 
work. 

If  at  least  one  of  the  functions  in  the  convolution  consists  entirely  of  exponentials, 
then  the  integration  can  be  performed  recursively.  This  can  result  in  significant  savings  of 
computer  time.  Several  authors  approximated  the  time-domain  impulse  responses  by  finite 
sums  of  exponentials  [25]-[26],  while  others  utilized  a  direct  frequency-domain  fitting  of 
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the  transfer  functions  by  rational  polynomials  [27] -[28].  With  direct  frequency-domain 
fitting,  the  need  to  perform  an  inverse  transform  is  eliminated. 

In  the  past  decade  there  have  been  a  number  of  attempts  to  solve  the  most  general 
problem  of  lossy  multiconductor  transmission  lines  with  nonlinear  terminations  [29]-[32]. 
Triezenberg  [29]  developed  a  state-variable  transmission  line  model  by  approximating  the 
spatial  variation  of  the  line  with  a  finite  Fourier  cosine  series.  Gruodis  and  Chang  [30] 
combined  the  decoupled-mode  transformation,  the  method  of  characteristics,  and  the  ratio¬ 
nal  polynomial  approximations  together  with  a  state  variable  solution.  Djordjevic  et  al.  [31] 
used  quasi-matched  Y-parameters  to  represent  the  multiconductor  lines.  The  introduction 
of  the  quasi-matched  parameters  necessitated  the  adding  of  negative  resistances  to  the 
terminal  networks.  Schutt-Aine  [32]  avoided  this  problem  by  using  scattering  parameters 
(S  parameters)  to  characterize  the  line.  Although  all  of  these  techniques  are  theoretically 
correct,  results  have  only  been  shown  for  the  simplest  case  of  two  symmetric  conductors 
over  a  ground  plane.  In  this  particular  case,  the  modal  transformation  matrices  result  in  an 
even  and  odd  mode  no  matter  what  the  frequency;  as  a  consequence,  the  most  general  case 
is  avoided. 

In  this  study  time-domain  models  for  lossy,  linear,  uniform  multiconductor  trans¬ 
mission  lines  are  developed  based  on  discretized  Thevenin  and  Norton  equivalent  circuits. 
The  analysis  is  quasi-TEM  in  that  the  models  are  derived  from  the  (possibly  frequency- 
dependent)  per-unit-length  resistance,  inductance,  conductance,  and  capacitance  matrices 
which  are  assumed  to  be  available.  These  models  are  then  incorporated  into  a  general  cir¬ 
cuit  analysis  program  for  simulating  dynamic  nonlinear  circuits.  Using  the  per-unit-length 
line  parameters,  the  exponential  propagation  matrix  and  the  characteristic  impedance  or 
admittance  matrix  are  calculated  for  a  range  of  frequencies.  The  exponential  propagation 
matrix  is  first  calculated  in  the  decoupled-mode  form  (along  with  the  modal  transformation 
matrices),  but  is  then  recombined  in  the  frequency  domain  to  the  full-matrix  form.  In  the 
full-matrix  form  of  the  exponential  propagation  matrix,  all  impulse  responses  have  a 
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physical  meaning.  The  fast  Fourier  transform  is  utilized  to  conven  the  transfer  functions  to 
time-domain  impulse  responses;  then,  direct  integration  of  the  convolution  integrals  is 
carried  out  during  the  time-stepping  procedure.  Results  are  shown  for  networks  of  suffi¬ 
cient  generality  to  fully  test  the  method. 

In  Chapter  2  a  general  method  for  computer-aided  analysis  of  dynamic  nonlinear 
circuits  is  presented.  The  theoretical  and  numerical  development  of  the  time-domain 
models  for  lossy  multiconductor  transmission  lines  is  described  in  Chapter  3.  In  order  to 
verify  the  simulation  technique,  a  computer  program  compatible  with  the  analysis  method 
of  Chapter  2  has  been  written  incorporating  the  lossy  line  models  of  Chapter  3.  This  pro¬ 
gram  is  demonstrated  in  Chapter  4  for  several  circuit  configurations.  Conclusions  drawn 
from  this  study  are  discussed  in  Chapter  5. 
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CHAPTER  2 

COMPUTER-AIDED  ANALYSIS 
OF  DYNAMIC  NONLINEAR  CIRCUITS 

2.1  Introduction 

Lumped-element  networks  are  completely  characterized  by  Kirchhoffs  current  law 
(KCL),  Kirchhoff’s  voltage  law  (KVL),  and  the  branch  element  equations.  In  general, 
when  energy-storage  and  nonlinear  elements  are  present  in  a  circuit  (i.e.,  a  dynamic  non¬ 
linear  circuit),  a  set  of  nonlinear  implicit  differential-algebraic  equations  results.  Unfortu¬ 
nately,  because  there  are  usually  no  closed-form  solutions  to  these  equations,  numerical 
methods  must  be  used. 

In  this  chapter  common  techniques  are  combined  to  form  a  method  of  computer- 
aided  analysis  of  dynamic  nonlinear  circuits.  The  objective  is:  if  circuits  that  contain  lossy 
multiconductor  transmission  lines  can  be  included  in  this  method,  then  they  can  be  included 
in  similar  circuit  analysis  programs  that  currently  exist.  The  technique  presented  herein 
utilizes  modified  nodal  analysis  to  formulate  the  circuit  equations,  the  trapezoidal  algorithm 
(with  a  fixed  time  step)  for  discrete  integration,  and  the  Newton-Raphson  algorithm  to  iter¬ 
atively  solve  systems  of  nonlinear  equations. 

2.2  Modified  Nodal  Analysis  of  Dynamic  Nonlinear  Circuits 

Modified  nodal  analysis  (MNA)  [33]  is  an  extension  of  nodal  analysis.  It  is  a 
completely  general  method  for  analyzing  dynamic  nonlinear  networks.  In  nodal  analysis 
KCL  is  formulated,  using  node  voltages  as  variables,  at  each  circuit  node  except  for  an 
arbitrarily  chosen  reference  node.  This  is  possible  only  if  the  current  through  each  element 
is  expressible  as  a  function  of  the  voltage  across  it.  Voltage  sources  and  nonlinear  current- 
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controlled  elements,  for  instance,  do  not  meet  this  requirement.  Therefore,  in  its  most 
basic  form,  nodal  analysis  may  be  inadequate  for  even  very  simple  circuits.  In  MNA, 
whenever  an  element  is  encountered  that  is  not  voltage  controlled,  the  current  through  the 
element  is  introduced  as  a  variable  in  the  node  equation,  and  the  branch  equation  for  that 
element  is  appended  to  the  system  of  equations  to  be  solved.  The  end  effect  is  a  set  of 
equations  in  which  the  unknowns  are  the  node  voltages  and  some  selected  branch  currents. 

For  circuits  that  contain  nonlinear  and  dynamic  elements,  the  network  equations 
have  the  following  functional  form  [34]: 

=  0  (2.1a) 

f2ix^,X2,...^m^\JC2,.--,Xi,t)  =  0  (2.1b) 


=  0  (2.1m) 

where  t  is  the  temporal  variable  and  the  x,’s  are  the  unknown  circuit  variables.  Also,  /  <  m 
and  the  dot  above  the  x,’s  indicates  the  derivative  with  respect  to  time.  Defining  vectors  x 
and  X  as 

\  =  [xx  X2  ■  ■  ■  xj^ 

X  =  [xi  X2  ■  ■  ■ 

where  the  T  indicates  transpose,  Eq.  (2. 1)  can  be  expressed  concisely  as 

f(x,x,t)  =  0  (2.2) 

in  which  0  is  a  vector  of  zeros. 

Chua  and  Lin  [34]  provide  detailed  procedures  for  solving  the  implicit  differential- 
algebraic  (IDA)  system  specified  by  Eq.  (2.2).  A  brief  overview  is  given  here.  Suppose 
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that  the  solution  x(0  of  Eq.  (2.2)  had  been  found  at  time  t  =  t^_  =  kh  (and  all  previous  times 
tj  for  j  =  0,l,2,...,j(:-l),  where  h  is  the  time  step.  Denoting  x{tj)  by  Xj,  the  solution  at 
t  =  must  satisfy 

f(Xi^i,x(ti+i),ri+i)  =  0  (2.3) 

The  trapezoidal  algorithm  can  be  used  to  approximate  the  next  value  of  x(4+i)»  die  time 
derivative  of  x(r)  evaluated  at  r  =  in  terms  of  xCt^^i),  x(ti^),  and  past  values  of  x(r). 
However,  this  integration  technique  cannot  be  applied  to  Eq.  (2.3)  in  its  present  form; 
instead,  the  trapezoidal  rule  must  be  applied  directly  to  the  branch  equations.  Introducing 
the  vector  u  of  branch  currents  and  voltages,  the  branch  equations  can  be  expressed  ana¬ 
lytically  as 

ii  =  g(u,r)  (2.4) 


where 

T 

U  =  [wi  •  •  •  M2/}] 

T 

il  =  [Mi  M2  •  •  •  M2/j] 

and  is  the  number  of  branches.  Applying  the  trapezoidal  integration  algorithm  to  the 
system  of  differential  equations  specified  by  Eq.  (2.4)  results  in  the  following  formula: 

Uyk+I  =  Uik  +  Ytg(Ufc,ffc)  +  g(ut+i,ri+i)l  (2.5) 

With  the  branch  equations  in  the  form  of  Eq.  (2.5),  the  modified  nodal  equations  are  recast 
with  x^^.1  as  the  unknown  variable.  This  results  in  x{ti^^.i)  being  replaced  by  x^^.j  and  x^, 
and  functions  that  depend  on  them.  Since  the  variable  x^^  is  known  at  time  Eq.  (2.3) 
may  then  be  rewritten  as 


(2.6) 
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Observe  that  Eq.  (2  6)  is  a  system  of  nonlinear  algebraic  equations  in  terms  of  the 
unknown  variable  The  most  common  technique  for  solving  systems  of  nonlinear 
equations  is  the  Newton-Raphson  (N-R)  algorithm.  Using  the  N-R  algorithm,  Xi^.j  is 
iterated  to  a  solution  via  the  following  recursive  relation: 

-  f(x2!i,r*+i)  (2.7) 


where  the  superscript  indicates  the  current  iteration  and 


9f(x,r) 

dx 


is  the  Jacobian  matrix  of  f(x,r)  evaluated  at  x  =  x^^j.  Equation  (2.7)  is  a  linear  algebraic 
system  of  equations  of  the  form  Ax  =  b,  where  A  is  a  constant  matrix,  b  is  a  constant 
vector,  and  x  is  a  vector  of  unknowns.  There  are  several  efficient  numerical  algorithms  for 
its  solution. 

Summarizing,  the  method  just  outlined  requires  three  steps.  First,  the  circuit  equa¬ 
tions  are  formulated  using  the  MNA  procedure.  This  results  in  a  nonlinear  implicit  differ¬ 
ential-algebraic  system  of  equations.  Next,  the  IDA  equations  are  discretized  by  applying 
the  trapezoidal  rule  to  the  branch  equations.  This  leads  to  a  system  of  nonlinear  algebraic 
equations  that  must  be  solved  at  each  time  step.  And  last,  the  nonlinear  equations  are  lin¬ 
earized  using  the  N-R  technique.  This  leaves  a  set  of  linear  algebraic  equations  to  be 
solved  at  each  iteration. 


2.3  Companion  Models  and  Element  Stamps 

In  computer  implementation  of  MNA,  KCL  and  the  branch  equations  are  not 
assembled  one  at  a  time  as  they  are  in  hand  calculations.  Instead  each  element  of  the  circuit 
is  modeled  by  first  discretizing  (if  necessary)  and  then  linearizing  (again  if  necessary)  the 
branch  equation(s)  for  that  element.  The  resulting  equivalent  circuit,  which  depends  on  the 


given  integration  and  nonlinear  iteration  schemes,  is  commonly  referred  to  as  the  element’s 
companion  model.  The  discretizing  of  the  elements  generates  a  nonlinear  static  (resistive) 
circuit,  instead  of  the  nonlinear  dynamic  circuit  originally  given.  The  linearization  changes 
the  nonlinear  resistive  circuit  into  a  linear  resistive  circuit.  The  network  equations  are  then 
formulated  using  the  rules  previously  outlined  for  MNA.  Since  the  network  now  consists 
solely  of  sources  and  resistors,  the  circuit  equations  are  of  the  form  Ax  =  b.  The  overall 
solution  of  the  network  amounts  to  solving  a  sequence  of  nonlinear  resistive  circuits  by 
iterating  a  sequence  of  linear  resistive  circuits.  The  derivation  of  companion  models  for  a 
number  of  common  circuit  elements  will  now  be  considered. 

The  branch  equations  for  the  capaci- 

:  i 

where  C  is  the  capacitance,  v  is  the  voltage  across  the  capacitor,  q  is  the  charge  on  one  of 
its  plates,  and  /  is  the  current  through  the  capacitor  from  node  y  to  node  z.  Using  the  trap¬ 
ezoidal  method  to  discretize  (7  =  /  at  time  results  in 

h  . 

+  h+\> 

Substituting  Cv  for  q  and  solving  for  yield 

~  ^  (2.8) 

As  shown  in  Fig.  2.1(b),  Eq.  (2.8)  can  be  realized  as  a  resistor  with  conductance  2Clh  in 
parallel  with  a  current  source  of  value  ICv/Jh  +  i/^.  The  conductance  value  is  Fixed  because 
the  time  step  h  is  constant,  but  due  to  the  memory  inherent  in  this  energy- storage  element. 
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Figure  2. 1  Linear  capacitor:  (a)  symbol  and  definitions,  (b)  companion  model  at  time 
tfc+i,  (c)  element  stamp,  (d)  element  stamp  with  branch  current  as  a  variable. 
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the  current  source  depends  on  the  previous  capacitor  voltage  and  current.  This  means  that 
the  current  source  will  have  to  be  updated  at  each  time  step. 

Following  the  rules  of  MNA,  one  can  visualize  writing  KCL  at  nodes  y  and  z  to 
arrive  at  the  nonzero  entries,  shown  in  Fig.  2.1(c),  of  the  coefficient  matrix  A  and  source 
vector  b.  Since  the  current  through  the  capacitor  is  expressible  as  a  function  of  the  voltage 
across  it,  it  was  not  necessary  to  introduce  the  current  as  a  variable  when  implementing 
KCL.  However,  it  is  sometimes  convenient  to  do  so,  particularly  when  the  current  is 
desired  as  an  output  variable,  or  when  it  is  needed  for  a  current-controlled  element  In  such 
cases  KCL  is  written  at  each  node  with  the  current  as  the  variable,  and  then  the  discretized 
branch  equation,  Eq.  (2.8),  is  appended  to  get  the  nonzero  entries  of  A  and  b  indicated  in 
Fig.  2.1(d). 

An  analogous  procedure  is  used  to  derive  the  companion  model  for  a  linear  inductor 
and  its  corresponding  element  stamp.  The  necessary  equations  for  the  inductor  depicted  in 
Fig.  2.2(a)  are 


<l>  =  Li 


d(p 

IF 


V 


where  L  is  the  inductance,  v  is  the  voltage  across  the  inductor,  (p  is  the  flux  linking  the 
inductor  coils,  and  i  is  the  current  through  the  inductor  from  node  y  to  node  z.  Following 
steps  identical  to  those  just  used  for  the  capacitor,  the  discretized  branch  equation  for  the 
inductor  at  time  is  easily  found  to  be 


(2.9) 


This  equation  leads  to  the  companion  model  illustrated  in  Fig.  2.2(b).  Since  a  voltage 
source  is  contained  in  the  model,  the  stamp  for  this  element  is  derived  by  writing  KCL  at 
each  terminal  of  the  element  using  the  current  through  the  inductor  as  the  variable,  and  then 
Eq.  (2.9)  is  added  to  get  the  A  and  b  entries  shown  in  Fig.  2.2(c). 
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2.3.2  Controlled  sources  and  coupled  elements 

Controlled  sources  are  easily  included  in  computer-aided  circuit  analysis  programs 
that  utilize  the  modified  nodal  approach.  The  current-controlled  current  source  (CCCS)  is 
employed  as  an  example  here. 

The  symbol  for  the  CCCS  is  illustrated  in  Fig.  2.3(a).  The  branch  equation  for  this 
element  is  /  =  a  i^,  where  is  the  controlling  current  through  some  element  a,  and  a  is  the 
gain  of  the  source.  The  CCCS  is  a  memoryless  element;  consequently,  its  companion 
model  at  any  time  step  is  just  the  source  itself.  This  model  is  shown  in  Fig.  2.3(b)  and  the 
corresponding  element  stamps  are  depicted  in  Figs.  2.3(c)  and  2.3(d).  When  constructing 
these  stamps,  the  controlling  current  was  assumed  to  be  available  as  a  circuit  variable. 
The  stamps  of  other  types  of  controlled  sources  are  derived  in  much  the  same  manner. 

Coupled  elements  can  generally  be  characterized  by  some  combination  of  linear 
resistors  and  controlled  sources.  Take  for  instance  the  coupled  inductors  illustrated  in  Fig. 
2.4(a).  The  equations  for  this  element  can  be  written  compactly  in  matrix  form  as 


where  v  is  a  vector  of  branch  voltages,  (p  a  vector  of  related  fluxes,  i  a  vector  of  branch 
currents,  and  L  the  inductance  matrix.  Discretizing  9  =  v  at  time  substituting  Li  for 
(p,  and  solving  for  yield 


'k+\ 


=  TLifc+i  -  -fLit+v 


(2.10) 


Expanding  Eq.  (2.10),  the  following  is  obtained: 


2  .  2  .  r2  .  2 

''i  k+\  =  ■^^11*1  k+\  +  "J^^hk+i  ~  k  +  i  k 

2  .  2  .  (2.2. 

^2  k-<-\  -  k+\  +  k+\  ~  I  k  +  '11^22^2  k'^  ^2k 


(2.11a) 

(2.11b) 
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Figure  2.3  Linear  current-controlled  current  source:  (a)  symbol  and  definitions,  (b)  com¬ 
panion  model  at  time  (c)  element  stamp. 
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Figure  2.4  (cont.)  Linear  coupled  inductors:  (c)  element  stamp. 
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From  Eq.  (2.1 1)  one  can  see  that  the  diagonal  elements  of  L  multiplied  by  Hh  are  resis¬ 
tance  values,  while  the  off-diagonal  elements  multiplied  by  Hh  are  transresistance  values  of 
gain  for  current-controlled  voltage  sources.  The  last  terms  on  the  right-hand  sides  of  Eqs. 
(2. 1  la)  and  (2.1  lb)  can  be  represented  by  independent  voltage  sources. 

The  companion  model  for  the  coupled  inductors  is  shown  in  Fig.  2.4(b).  The 
stamp  is  constructed  in  a  manner  similar  to  that  for  a  single  inductor  and  is  presented  in 
Fig.  2.4(c). 

2.3.3  The  nonlinear  capacitor 

Nonlinear  resistors  and  capacitors  are  the  nonlinear  elements  most  needed  when 
modeling  digital  devices.  Since  the  nonlinear  capacitor  actually  contains  a  nonlinear  resis¬ 
tor  in  its  companion  model,  the  capacitor  will  be  used  as  the  demonstration  element  here. 

Figure  2.5(a)  shows  a  nonlinear  capacitor.  Most  nonlinear  capacitors  encountered 
are  voltage  controlled  and  can  be  represented  analytically  by  the  following  equations: 


where  q,  v,  and  i  are  charge,  voltage,  and  current,  respectively,  and  q{v)  indicates  that  the 
charge  is  some  nonlinear  function  of  the  voltage.  Discretizing  the  time  derivative  of  the 
charge  using  the  trapezoidal  rule  and  solving  for  result  in 

‘fc+i  =  ^  ^2.12) 

This  equation  can  be  represented  by  the  equivalent  circuit  shown  in  Fig.  2.5(b).  It  is  simi¬ 
lar  to  the  model  for  the  linear  capacitor  of  Fig,  2. 1(b),  except  that  now  the  resistor  is  non¬ 
linear. 

From  Eq.  (2.12)  and  Fig.  2.5(b),  one  can  see  that  the  current  through  the  nonlinear 
resistor  obeys  the  equation  i/j  =  2q(v)/h.  A  sample  v-ip  curve  for  this  function  has  been 
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drawn  in  Fig.  2.6(a).  Suppose  that  the  solution  for  the  network  had  been  found  at  time  tj^. 
This  result  can  be  used  to  initiate  an  iteration  scheme  to  find  the  solution  at  time  After 
the  (/)th  iteration,  the  solution  for  the  nonlinear  resistor  will  be  a  particular  point  on  its  v-/)j 
curve.  Designating  this  point  as  on  the  sample  curve  of  Fig.  2.6(a),  the  operation  of 
the  resistor  can  be  approximated  by  the  line  drawn  tangent  to  the  curve  at  that  point.  The 
equation  for  this  line  is 

(2.13) 


where 


2  Bqjv) 
h  3v 


is  the  slope  at  point  and  is  the  i/j  intercept.  Solving  Eq.  (2.13)  for  this  intercept 
yields 

!  0  )  _  ;  0  )  _ 

‘S  k+l  -  k+\  Sk-¥\'^k+\ 

Now  substituting  back  into  Eq.  (2.13),  and  updating  the  iteration  count  on  the 
unknown  variables  iRk,+\  and  the  following  recursive  formula  is  obtained: 


;  0+1)  _  „0)  .0+1)  ,  r 

1/?  fe+1  -  8k+l^k+l  +  ‘ 


0) 

s  ,t+i 


(2.14) 


Equation  (2.14)  can  be  modeled  by  the  parallel  combination  of  resistor  and  current  source 
shown  in  Fig.  2.6(b).  Substituting  Eq.  (2.14)  into  Eq.  (2.12)  yields  (since  =  2q{v)/h) 

h^Li  ~  (2.15) 


Equation  (2.15)  leads  to  the  overall  companion  model  of  Fig.  2.7(a)  for  the  nonlinear 
capacitor.  This  model  characterizes  the  capacitor  for  iteration  (/+1)  at  time  The 
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Figure  2.7  Nonlinear  capacitor;  (a)  companion  model  for  iteration  (j+\)  at  time 
(b)  element  stamp,  (c)  element  stamp  with  branch  current  as  a  variable. 
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linearization  method  presented  here  can  be  shown  to  be  identical  to  the  Newton-Raphson 
method  discussed  in  Section  2.2  [34]. 

Comparing  the  companion  model  for  the  nonlinear  capacitor  here  with  that  of  the 
linear  capacitor  of  Fig.  2.1(b),  it  can  be  seen  that  they  have  the  same  form.  This  being  the 
case,  the  element  stamps  for  the  nonlinear  capacitor  are  easily  constructed  and  are  illustrated 
in  Figs.  2.7(b)  and  2.7(c).  When  using  the  stamp  of  Fig.  2.7(b),  convergence  of  the  volt¬ 
age  within  the  prescribed  limits  does  not  always  mean  that  sufficient  convergence  of  the 
current  has  occurred.  For  this  reason,  the  stamp  of  Fig.  2.7(c)  is  generally  preferred  for 
nonlinear  elements  since  the  current  is  available  as  a  variable  and  may  be  subjected  to  the 
same  convergence  requirements  as  the  voltage. 

2.4  Conclusions 

In  this  chapter,  a  commonly  used  and  completely  general  method  of  computer-aided 
circuit  analysis  was  presented.  Implementation  of  a  particular  element  in  a  circuit  analysis 
program  consists  of  developing  a  companion  model  for  that  element  and  substituting  the 
element’s  stamp  into  the  coefficient  matrix  A  and  source  vector  b.  The  overall  solution  of  a 
network  amounts  to  solving  a  sequence  of  nonlinear  resistive  circuits  by  iterating  a  se¬ 
quence  of  linear  resistive  circuits.  In  the  following  chapter,  two  different  companion 
models  and  their  corresponding  element  stamps  are  developed  for  a  lossy  multiconductor 
transmission  line  in  order  that  it  may  be  included  in  a  circuit  analysis  program  of  the  type 
discussed  in  this  chapter. 
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CHAPTER  3 

LOSSY  MULTICONDUCTOR  TRANSMISSION  LINES 

3.1  Introduction 

The  analysis  of  nonlinear  circuits,  in  the  general  case,  can  be  performed  only  in  the 
time  domain.  On  the  other  hand,  the  analysis  of  lossy  transmission  lines,  particularly  those 
with  frequency-dependent  parameters,  is  handled  most  conveniently  in  the  frequency 
domain.  Therefore,  in  order  to  devise  a  method  of  analysis  for  lossy  multiconductor 
transmission  lines  in  nonlinear  circuits,  it  is  necessary  to  combine  solutions  from  the  two 
domains.  Assuming  that  the  transmission  lines  are  linear,  the  inverse  Fourier  transform 
can  be  used  to  convert  frequency-domain  characterizations  of  the  lines  to  the  time  domain. 
Boundary  conditions  at  the  interface  of  the  lossy  distributed  lines  and  lumped  nonlinear 
networks  can  then  be  matched  at  discrete  time  steps.  [31  ] 

In  this  chapter,  the  lossy  multiconductor  transmission  line  equations  are  first  solved 
in  the  frequency  domain;  then,  the  inverse  Fourier  transform  is  utilized  to  switch  to  the  time 
domain  where  two  different  companion  models  are  developed  to  represent  the  multicon¬ 
ductor  line  at  its  endpoints.  One  of  the  models  consists  of  a  discretized  Thevenin  equiva¬ 
lent  circuit  for  each  end  of  the  line,  while  the  other  model  consists  of  discretized  Norton 
equivalent  circuits.  Depending  on  the  properties  of  the  line,  it  may  be  advantageous  to  use 
one  or  the  other.  In  addition,  numerical  issues  are  discussed  for  calculating  and  imple¬ 
menting  the  models  in  a  general  circuit  analysis  program. 

3.2  Frequency-Domain  Solutions 

A  lossy  (n+l)-conductor  transmission  ''np  of  length  'i  is  shr\wn  'n  Fig,  3.1.  The 
telegrapher’s  equations  for  this  transmission  line  can  be  expressed  in  matrix  form  as  [35] 
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av(x,r)  ^  ^ 


dx 

3I(££) 


=  GV(jc./)  +  C 


avcc,f) 

ar 


(3.1a) 

(3.1b) 


where  x  and  t  are  the  spacial  and  temporal  variables,  respectively;  V(x,r)  is  an  n  x  1  vector 
of  line  voltages;  I(x,r)  is  an  nx  1  vector  of  line  currents;  and  R,  L,  G,  and  C  are  the  n  x  n 
per-unit-length  resistance,  inductance,  conductance,  and  capacitance  matrices,  respectively. 
The  yth  element  of  V(x,r),  Vj(x,t),  and  the  yth  element  of  I(x,r),  Ij{x,t),  are  the  voltage  and 
current,  respectively,  on  line  j.  In  the  formulation  above,  the  line-parameter  matrices  are 
constant;  this  is  of  limited  usefulness  since  for  even  relatively  slowly  varying  signals  the 
relationship  between  voltage  and  current  on  most  transmission  systems  is  actually  more 
complicated  [35].  Once  Eq.  (3.1)  has  been  transformed  into  the  frequency  domain,  this 
restriction  on  the  time  rate  of  change  of  voltages  and  currents  can  be  overcome  in  many 
cases  by  allowing  the  line-parameter  matrices  R,  L,  G,  and  C  to  be  frequency  dependent. 
Proceeding  in  this  manner,  the  Fourier  transform  of  Eq.  (3.1)  with  respect  to  time  t  yields 

- —  =  (R  -I-  j(oL)l{x,co) 

- ^ —  =  (G  +  jcoC)y(x.CD) 

where  co  is  the  radian  frequency;  then,  explicitly  showing  the  frequency  dependence  of  the 
line-parameter  matrices  leads  to 

=  (R(<y)  +Ja}L(ct))]l(x,co)  =  Z(ty)I(x,<y)  (3.2a) 

=  [G(tu)  -I-  y£wC(Q))]V0c,£y)  =  Y(£u)V(;c,tu)  (3.2b) 

where  the  impedance  per  unit  length,  Z(<y)  =  R(t»)  +  Jq}L(co),  and  the  admittance  per  unit 
length,  Y(cu)  =  G(co)  +  jcoCico),  have  been  introduced.  As  a  consequence  of  reciprocity. 
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R((U),  L(ct)),  Gico),  and  C(co)  are  symmetric;  therefore,  so  are  Z(6))  and  Y(£o).  A  physical 
interpretation  of  R(tu),  L(<u),  G(0)),  and  Cioo)  is  given  later  in  Sec.  3.2.1.  The  calculation 
of  these  matrices  is  not  a  part  of  this  work.  In  this  study,  for  any  simulation,  values  for 
R(to),  L(tu),  Gico),  and  C(6))  will  be  assumed  to  be  available  for  each  transmission  line 
encountered;  they  can  be  measured  or  calculated  by  one  of  the  several  available  techniques 
[36]-[41]. 

In  order  to  characterize  muldconductor  systems  by  the  line-parameter  matrices,  the 
modes  of  propagation  must  be  transverse  electromagnetic  (TEM)  or  quasi-TEM.  In  multi¬ 
conductor  transmission  line  systems  with  conductor  losses,  whether  the  systems  are  con¬ 
structed  with  homogeneous  or  nonhomogeneous  dielectrics,  the  existence  of  purely  TEM 
modes  is  nut  possible;  however,  as  long  as  the  losses  are  small,  quasi-TEM  modes  can 
propagate.  In  fact,  as  long  as  the  wavelength  of  the  highest  frequency  under  consideration 
is  long  compared  to  the  greatest  cross-sectional  distance  of  the  system,  only  quasi-TEM 
modes  can  propagate.  In  addition  to  these  implicit  constraints,  it  also  will  be  assumed  that 
the  line  is  of  uniform  cross  section  throughout  its  length;  that  is,  R,  L,  G,  and  C  are 
assumed  to  be  constant  with  respect  to  the  spacial  variable  x. 

Proceeding  now  to  the  solution  of  Eq.  (3.2),  substituting  Eq.  (3.2b)  into  the 
derivative  of  Eq.  (3.2a)  and  Eq.  (3.2a)  into  the  derivative  of  Eq.  (3.2b)  yields  one  equation 
in  terms  of  the  line  voltages  and  another  in  terms  of  the  line  currents: 

—  =  ZYV  =  r^V  (3.3a) 

— T  =  YZI  =  Y^Z^I  =  (ZY)^I  =  (r^)^I  (3.3b) 

a^^ 

where  the  propagation  matrix  f  has  been  introduced.  In  Eq.  (3.3a),  =  ZY  has  been 

detined,  and  due  to  the  symmetry  of  Z  and  Y,  YZ  then  equals  (F")  as  shown  in  Eq. 
(3.3b).  In  general,  even  though  Z  and  Y  are  symmetric,  F^  is  not.  If  the  existence  of  two 
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nonsingular  square  matrices  E  and  J  is  assumed,  then  left-multiplying  both  sides  of  Eq. 
(3.3a)  by  E"^  and  both  sides  of  Eq.  (3.3b)  by  J"'  results  in 

^  ^  =  E'V^EE'^V  (3.4a) 

^^  =  rkryjJ-4  (3.4b) 


where  the  identity  matrix  has  been  inserted  in  the  right-hand  sides  of  Eqs.  (3.4a)  and  (3.4b) 
as  EE~^  and  JJ~\  respectively.  The  matrix  products  E“'r^E  and  J~\r^)^J  in  Eq.  (3.4) 
represent  similarity  transformations  of  and  (F^)^,  respectively.  If  there  exists  n  inde- 
pendent  eigenvectors  of  F  and  n  independent  eigenvectors  of  (F  )  ,  and  the  columns  of  E 
and  J  are  chosen,  respectively,  to  be  these  eigenvectors,  then  the  aforementioned  similarity 
transformations  result  in  diagonal  matrices.  The  existence  of  eigenvectors  that  diagonalize 
1  and  (F  )  will  be  assumed  throughout  this  investigation.  In  fact,  it  is  well  known  that  a 
matrix  and  its  transpose  are  similar  to  the  same  diagonal  matrix  [42];  designating  this  matrix 
as  A""  yields 

A^  =  E"’F^E  =  rkFyj  (3.5) 

The  diagonal  elements  of  A  ate  the  eigenvalues  of  F  and  (P)  .  Their  square  roots  are  the 
complex  quasi-TEM  propagation  modes  of  the  lossy  multiconductor  transmission  line;  for 
this  reason,  A  will  be  referred  to  as  the  modal  propagation  matrix.  The  yth  columns  of  E 
and  J  are  the  voltage  and  current  eigenvectors,  respectively,  associated  with  the  yth  eigen- 
value  of  A  .  Now  substituting  Eq.  (3.5)  into  Eq.  (3.4)  and  defining  the  modal  variables 
=  E~W  and  lead  to 


dx^ 


(3,6a) 
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(3.6b) 


The  modal  voltage  equation,  Eq.  (3.6a),  and  the  modal  current  equation,  Eq.  (3.6b),  can  be 
solved  independently  of  one  another.  First,  the  solution  to  the  voltage  equation  will  be 
derived.  Since  A  is  diagonal,  the  individual  equations  in  Eq.  (3.6a)  are  decoupled.  This 
means,  for  instance,  that  modal  voltage  the  yth  element  of  V„,  depends  only  on  the 
derivative  of  Vj^  and  no  other  modal  voltages;  therefore,  Eq.  (3.6a)  is  actually  n  separate 
scalar  equations  of  the  form 


where  A^  is  the  y'th  diagonal  element  of  A  and  is  the  propagation  constant  for  mode  j.  The 
solution  of  this  equation  is 

V-  =  +  eVy  - 

'  jm  *  jm  ^  '  jm 

where  and  Vj^~  are  scalar  constants  that  depend  on  the  boundary  conditions.  Conse¬ 
quently,  the  solution  of  Eq.  (3.6a)  can  be  written  in  matrix  form  as 

=  +  (3.7) 

where  e~^,  the  modal  exponential  propagation  matrix,  is  a  diagonal  matrix  whose  yth 
diagonal  element  is  and  and  are  vector  constants.  Multiplying  each  side  of 
Eq.  (3.7)  on  the  left  by  E  and  inserting  the  identity  matrix  as  E~'e  in  each  term  on  the 
right-hand  side  result  in 


EV^  =  Ee''^E"‘EV,;  +  Ee^E~’EV  “ 
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Making  use  of  a  previous  definition  relating  modal  voltages  and  line  voltages  [V  =  EV^] 
and  the  fact  from  matrix  algebra  that  Ee~^E~^  =  *  [42]  yields 

where  line-voltage  constants  V'*’  =  EV^"^  and  V~  =  EV^“  have  been  defined.  Now  utiliz¬ 
ing  the  definition  of  A  from  Eq.  (3.5),  the  solution  for  the  voltage  as  a  function  of  x  is 
obtained: 

V  =  -H  (3.8) 

where  e~^^,  the  exponential  propagation  matrix,  has  been  introduced.  The  first  term  on  the 
right-hand  side  of  Eq.  (3.8)  describes  a  vector  of  positive  traveling  waves,  while  the 
second  term  represents  a  vector  of  negative  traveling  waves.  To  find  the  solution  for  the 
current,  Eq.  (3.8)  is  substituted  into  Eq.  (3.2a).  Taking  the  derivative  and  rearranging 
some  terms  lead  to 

ZqI  =  (3.9) 

where 

Zo  =  r"'z  (3.10) 

is  called  the  characteristic  impedance  matrix.  From  Eqs.  (3.8)  and  (3.9),  it  can  be  seen  that 
Zq  relates  the  positive  and  negative  traveling  voltage  waves  to  the  positive  and  negative 
traveling  current  waves,  respectively.  As  shown  in  Sec.  3.2.2,  the  characteristic 
impedance  matrix  Zq  is  symmetric.  Equations  (3.8)  and  (3.9)  are  the  solution  to  the  lossy 
multiconductor  iransmission  line  equation,  Eq.  (3.2).  Note  that  the  solution  has  2n  inde¬ 
pendent  constants  V"*"  and  V~  that  depend  on  the  boundary  conditions. 
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A  different  form  of  the  solution  of  Eq.  (3.2)  can  be  obtained  by  solving  the  modal 
current  equation,  Eq.  (3.6b).  If  this  is  done,  the  solutions  for  the  current  and  voltage  are 
obtained  as 

I  =  e"‘'V-/V  (3.11a) 

YoV  =  e‘^V  +  e^V  (3.11b) 

where 

Yo  =  (rV'Y  (3.12) 

is  the  characteristic  admittance  matrix  and  and  I”  are  constant  vectors.  The  minus  sign 
preceding  the  second  term  in  Eq.  (3.1  la)  arises  from  the  fact  that  V  is  defined  as  positive 
in  the  opposite  direction  of  I^,  which  has  the  same  reference  direction  as  the  total  current  I 
(see  Fig.  3.1).  As  might  be  expected,  Yq  is  the  inverse  of  Zq.  This  is  shown  to  be  the 
case  in  Sec.  3.2.2.  Note  that  this  means  that  Yq  is  also  symmetric. 

Equations  (3.8),  (3.9),  (3.11a),  and  (3.11b)  will  be  used  in  Sec.  3.3  to  derive  two 
time-domain  models  for  a  lossy  multiconductor  transmission  line.  But  first,  in  the  next  two 
subsections,  a  discussion  of  the  line -parameter  matrices  and  several  auxiliary  proofs  are 
given. 


3.2.1  Interpretation  of  the  line-parameter  matrices 

A  physical  interpretation  of  R(<w),  L((y),  G(£y),  and  C(co)  may  be  obtained  by  gen¬ 
eralizing  the  definitions  of  their  corresponding  scalar-valued  parameters.  Referring  to  Fig. 
3.1,  say  that  line  j  has  been  charged  to  one  volt  (in  the  time-harmonic  sense)  with  all  other 
lines  held  at  a  voltage  of  zero;  then,  Cy,  the  ijih  entry  of  C,  is  equal  to  the  charge  per  unit 
length  induced  on  line  i.  The  off-diagonal  elements  of  C,  sometimes  called  electrostatic 
inductance  coefficients,  are  always  negative.  The  diagonal  entries,  on  the  other  hand,  are 
always  positive  and  are  the  self-capacitances  of  the  individual  lines  of  the  n-line  system. 
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The  conductance  matrix  G  can  be  analyzed  in  a  similar  manner.  Again,  with  line  j  charged 
to  one  volt  and  all  other  lines  at  zero  volts,  the  net  current  leaving  line  i  per  unit  length  is 
entry  of  G.  As  is  the  case  with  the  capacitance  entries,  Gy  <  0  for  i  ^  j,  while  Gjj  >  0. 
As  one  might  expect,  an  entry  in  L,  such  as  L,y,  measures  the  flux  per  unit  length  linking 
the  loop  formed  by  line  i  and  the  reference  conductor  due  to  a  one-amp  current  flowing  in 
line  j.  Of  course,  diagonal  entries  represent  self-inductances  while  off-diagonal  entries  are 
referred  to  as  mutual  inductances.  Entry  Rij  of  R  measures  the  voltage  drop  along  line  i  per 
unit  length  due  to  a  one-amp  current  flowing  in  line  J.  For  i  -  j,  this  is  simply  the  usual  line 
resistance.  For  i  j,  it  may  seem  impossible  that  current  in  a  conductor  can  cause  a  voltage 
drop  in  a  conductor  other  than  itself.  Strictly  speaking  it  cannot;  however,  just  as  the  skin 
effect  of  an  isolated  conductor  can  cause  a  redistribution  of  current  that  results  in  an  effec¬ 
tive  increase  in  resistance  of  that  conductor,  the  time-varying  current  in  a  line  can  also  cause 
a  redistribution  of  current  in  a  neighboring  line.  Therefore,  the  effective  resistance  of  a 
conductor  in  the  neighborhood  of  other  current-carrying  conductors  depends  on  those  cur¬ 
rents.  This  dependency  is  modeled  by  the  off-diagonal  elements  of  R. 

3.2,2  Proof  that  Zq  =  and  Yg  =  Zg“* 

That  Zg  is  symmetric  and  its  inverse  equal  to  Yg  are  commonly  used  facts;  never¬ 
theless,  complete  and  efficient  proofs  are  hard  to  find.  For  instance,  the  symmetry  of  Zg 
has  been  argued  based  on  reciprocity  [4]  even  though  such  reasoning  should  be  reserved 
for  fundamental  physical  relations  such  as  the  line-parameter  matrices;  other  proofs  have 
been  either  incomplete  [5]  or,  in  the  corresponding  lossless-line  case,  inefficient  [18]. 

First,  it  must  be  shown  that  the  matrix  product  E“^ZJ  commutes  with  A“\  that  is, 

E"'ZJA~'  =  A'^E'^ZJ 
2 

The  starting  point  is  the  definition  of  A  from  Eq.  (3.5): 

E"'r^E  =  J~'(r^)^J  =  A‘ 


(3.13) 
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Substituting  ZY  for  and  remembering  that  Z  and  Y  are  symmetric  yield 

E'^ZYE  =  rWzj  = 

"■1  2  —1 

Inserting  the  identity  matrix  as  JJ  into  one  definition  of  A  and  as  EE  into  the  other 
definition  results  in 

E'^ZJJ'^YE  =  J-^YEE'^ZJ  =  A^ 

Defining  Az  =  E'^ZJ  and  Ay  =  J  Ve,  the  preceding  equation  can  be  rewritten  as 

AzAy  =  AyAz  =  A^ 

Now,  since  Az  and  Ay  commute  with  each  other,  they  also  commute  with  their  product  A^ 

2  _2  1 
and  its  inverse  A~  .  Also,  since  A  is  diagonal,  its  square  root  A~  is  diagonal. 

Therefore,  Az  and  Ay  also  commute  with  A"\  and  Eq.  (3.13)  is  proved. 

To  show  that  Zq  =  Zq\  one  begins  with  Zq  =  F'^Z  [Eq.  (3.10)];  then, 

zj'  =  (r^z)^  =  z(rV‘  *  zja'^j'^  =  ee"^zja~’j"^ 

Using  Eq.  (3.13)  yields  the  desired  result: 

Zo^  =  EA'^E-'ZJJ'^  =  EA"^E~‘Z  =  r"'z  =  Zq 
Next,  multiplying  Zq  =  F'^Z  by  Yq  =  (F^^Y  [Eq.  (3.12)],  one  obtains 
ZqYo  =  r  ^Z(fV‘y  =  EA‘^E''ZJA"'J"^Y 
Using  Eq.  (3.13)  again  leads  to 

ZoYo  =  ea"'a~‘e~*zjj~*y  =  =  1 

where  1  is  the  identity  matrix,  and  so  Yq  =  Zq~^. 
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3.3  Time-Domain  Models 

In  this  section,  two  time-domain  models  for  a  lossy  multiconductor  transmission 
line  are  developed.  One  of  these  is  based  on  Eqs.  (3.8)  and  (3.9),  the  solution  to  the  modal 
voltage  equation,  and  results  in  a  discretized  Thevenin  equivalent  circuit  for  each  end  of  the 
transmission  line;  the  other  is  based  on  Eqs.  (3. 11a)  and  (3. 1 1  b),  the  solution  to  the  modal 
current  equation,  and  leads  to  discretized  Norton  equivalent  circuits.  The  Thevenin  equiv¬ 
alent  model  will  be  developed  first  Both  the  space  dependence  and  frequency/time  depen¬ 
dence  of  the  equations  will  be  needed;  consequently,  they  will  be  shown  explicitly  in  the 
following  derivation.  Equations  (3.8)  and  (3.9)  are  rewritten  here  for  convenience: 

V{x,03)  =  -H  e^^‘“^'V"(a))  (3. 14a) 

Zo(£o)ICc,tu)  =  e"''^‘“^^V^(a))  -  e^^‘"^*V(G))  (3. 14b) 

Defining  \\x,Q))  =  e"^^‘“^V^(<u)  and  V“(jc,<u)  =  e^^“^V~(£t))  and  evaluating  these  equa¬ 
tions  at  jc  =  0  yield 

V^(x,ct>)  =  (3. 15a) 

V-(x,cu)  =  e’'^"^^V'(0,tu)  (3.15b) 

Substituting  Eq.  (3.15)  into  Eq.  (3.14)  and  evaluating  the  result  atx  =  0  lead  to 

V(0,to)  =  V"(0,tu)  -I-  ^(O.ty)  (3.16a) 

Zo((y)I(0,ty)  =  V^(0,ty)  -  V(0,ty)  (3. 16b) 

From  Eq.(3.15b),  one  notices  that  V~(0,G))  can  be  expressed  as  where  d  is 

the  length  of  the  multiconductor  line.  If  this  expression  is  substituted  into  Eqs.  (3.16a)  and 
(3.16b),  then  adding  the  respective  sides  of  the  two  resulting  equations  and  rearranging 
terms  yield 
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V(0,G))  =  V(4/,fi))  +  Zo(a))I(0,£o)  (3. 17) 

Equation  (3.17)  expresses  the  near-end  (x  =  0)  voltage  in  terms  of  the  near-end  current,  the 
far-end  (x  =  d)  reflected  voltages,  and  the  transfer  function  matrices  Zo(£i>)  and  e 
Similarly,  starting  with  Eq.  (3.14)  evaluated  at  the  far  end: 

W(d,Q))  =  V(0,tt))  -f-  e^^“^^V(0,m)  (3. 18a) 

-Zo(6))I(d,6))  =  V(0,tu)  -  /^“^^V(0,ty)  (3.1 8b) 

an  equation  is  obtained  that  expresses  the  voltage  at  the  far  end  in  terms  of  the  far-end  cur¬ 
rents,  the  near-end  reflected  voltages,  and  the  transfer  function  matrices; 

\id,Q})  =  2e"^^“^V(0,ty)  -h  Zo(ty)I(d,£u)  (3.19) 

Note  that  a  negative  sign  appears  on  the  left-hand  side  of  Eq.  (3.18b)  where  its  seems  there 
should  be  none.  In  the  original  lossy  multiconductor  transmission  line  equations,  Eq. 
(3.2),  the  reference  direction  for  the  current  was  chosen  in  the  positive  x  direction  as 
depicted  in  Fig.  3.1.  For  symmetry  purposes,  it  is  more  convenient  to  choose  the  reference 
direction  for  the  current  at  each  end  of  the  line  to  be  into  the  line.  This  necessitates  the 
introduction  of  the  negative  sign  in  Eq.  (3.18b). 

Before  transforming  Eqs.  (3.17)  and  (3.19)  to  the  time  domain,  a  modification  will 
be  made  to  the  characteristic  impedance  matrix  Zq.  In  general,  as  0)  approaches  infinity, 
Zo(£u)  approaches  a  real  constant  matrix;  therefore,  Zq  will  be  separated  into  the  sum  of 
two  terms  as  shown  below: 

Zg(ty)  =  Zqj{co)  +  Rg  (3.20) 

where  Rg  is  defined  as  the  limit  of  Zg( cu)  as  co  goes  to  infinity: 

Rg  =  lim  Zg( tu)  (3.21) 

Ct>  ><» 
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Equations  (3.20)  and  (3.21)  can  be  taken  to  define  the  “shifted”  characteristic  impedance 
matrix  Zq^(£0).  The  time-domain  Thevenin  companion  model  for  the  lossy  multiconductor 
transmission  line  can  now  be  constructed.  Since  Eqs.  (3.17)  and  (3.19)  have  the  same 
form,  it  will  be  necessary  to  show  a  detailed  derivation  for  only  one  end  of  the  line. 
Working  with  the  far  end,  Eq.  (3.20)  is  substituted  into  Eq.  (3.19)  to  obtain 

\{d,Q))  =  -I-  Zo,ico)I(d,Q))  +  RoKd.ta)  (3.22) 

Switching  to  the  time  domain,  the  inverse  Fourier  transform  of  this  equation  is  taken  with 
respect  to  co  to  arrive  at 

y(d,[)  =  2H(d,t)  *  v^(0,r)  +  Zo,(/)  *  lid,t)  +  RqK^,:)  (3.23) 

where  the  asterisk  denotes  convolution  and 

in  which  9  denotes  the  Fourier  transform  with  respect  to  t.  Explicitly  showing  the  convo¬ 
lution  integral  evaluated  at  time  r  =  l^^ds  to 

\id,t„)  =  2  /*  ‘H(d,r)V^(0,rfc-T)dT  -t-  f  "Zo,{r)l{d,t,,-T)dr  +  RoKd.r^)  (3.24) 
Jo  Jo 

In  the  computer  implementation  of  this  equation,  the  integrals  must  be  approximated  by 
finite  sums.  Rewriting  Eq.  (3.24)  using  finite  sums  results  in 

k-\  ^1 

V(d,tk)  =  2^H(c/,q)V^(0,fjt-r;+i)Ar  +  ^Zo,(r,)I(tf,rfc-r/+i)Ar  RoUdj^)  (3.25) 
1=0  1=0 

The  vector  V^(0,r/),  shown  in  the  first  term  on  the  right-hand  side  of  the  above 

equation  consists  of  all  of  the  reflected  voltages  that  have  been  placed  on  the  line  at  the  near 
end  since  the  beginning  of  the  simulation.  In  Eq.  (3.24)  H(d,z)  and  Zoj(r)  are  impulse 
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response  matrices;  this  means,  for  instance,  that  H(rf,r;),  /=0,...,^-l,  in  Eq.  (3.25)  is  a 
pulse  response  matrix.  This  pulse  response  matrix  takes  the  pulses  V'^(0,r/)Ar,  /=!,.. .,k, 
that  are  placed  on  the  line  at  a:  =  0  and  propagates  them  down  the  line  to  the  far  end;  there¬ 
fore,  the  summation  in  the  first  term  is  just  the  positive-going  voltage  wave  vector  arriving 
at  X  =  d  at  time  that  is,  The  second  summation  on  the  right-hand  side  of  Eq. 

(3.25)  involves  the  product  of  the  pulse  response  matrix  Zoj(t/),  and  a  vector 

of  past  current  pulses  at  the  far  end  of  the  line.  It  results  in  a  “past-history”  voltage  source 
that  represents,  to  a  certain  extent,  the  history  at  the  boundary  between  the  line  and  its 
termination  at  the  far  end  since  the  beginning  of  the  simulation.  This  voltage  source  will  be 
designated  V;,(d,r^.  Summarizing,  the  following  variables  have  been  defined: 

fe-i 

V^(d,r^)  =  ^H(d,r/)V^(0,ri-r,+i)Ar  (3.26a) 

1=0 

k-l 

=  X  (3.26b) 

1=0 


Substituting  these  variables  into  Eq.  (3.25)  leads  to 

V(d,;^)  =  2V^(d,4)  +  V;,(d,/P  +  RoI(d,/i)  (3.27) 

This  equation  expresses  the  voltage  at  the  far  end  at  time  in  terms  of  the  far-end  currents, 
the  constant  resistance  matrix  Rq,  and  two  known  voltages;  it  is  the  time  domain  counter¬ 
part  of  Eq.  (3.19).  A  similar  derivation  leads  to  an  equation  for  the  near  end  of  the  line. 
This  equation  may  be  written  as 

V(0,rt)  =  2V(0,rP  +  ¥^(0,/^)  +  RqKO,:*)  (3.28) 


where 

k-\ 

V~(0,4)  =  ^H(d,r,)V“(d,r^-r,^i)Ar 
1=0 


(3.29a) 
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k-l 

VftCO.fi)  =  51z05(^i)K0,ri-r,+i)Af  (3.29b) 

1=0 

A  lumped-element  companion  model  for  the  lossy  multiconductor  line  can  now  be 
developed  from  Eqs.  (3.27)  and  (3.28).  The  first  two  terms  of  the  right-hand  side  of  each 
equation  can  be  represented  by  a  voltage  source,  while  the  third  term  indicates  that  the  cur¬ 
rents  at  each  end  of  the  line  pass  through  a  system  of  coupled  resistances.  The  value  of  the 
voltage  source  is  the  open-circuit  voltage  at  time  that  results  from  an  incident  wave.  The 
set  of  coupled  resistances  represented  by  the  matrix  Rq  is  the  instantaneous  resistance 
looking  back  into  an  infinite  lossy  multiconductor  transmission  line.  This  justifies  the  use 
of  the  term  Thevenin  companion  model  for  the  equivalent  circuit  illustrated  in  Fig.  3.2.  In 
this  figure,  the  resistance  matrices  are  shown  in  block-diagram  form.  These  matrices  can 
be  implemented  in  a  manner  identical  to  the  resistance  matrix  encountered  in  Chapter  2  for 
the  discretized  coupled  inductors:  the  diagonal  elements  of  the  matrix  can  be  modeled  as 
two-terminal  resistors,  and  the  off-diagonal  elements  as  current-controlled  voltage  sources. 
Since  the  model  for  the  lossy  multiconductor  transmission  line  can  be  expressed  as  an 
assemblage  of  lumped  elements,  it  can  be  easily  included  in  any  circuit  analysis  program. 
This  has  been  a  major  goal  of  this  research. 

The  element  stamp  for  the  far  end  of  the  line  is  shown  in  Fig.  3.3.  Since  the  model 
contains  voltage  sources  in  each  line,  the  line  currents  are  also  introduced  as  variables  via 
the  branch  equations,  Eq.  (3.27).  These  currents  are  needed  in  any  case  because  their 
values  must  be  saved  so  that  they  are  available  for  convolution  with  Zo/r).  The  topology 
of  the  stamp  for  the  near  end  of  the  line  is  identical  to  that  for  the  far  end  and,  therefore,  is 
not  shown. 

A  Norton  companion  model  for  the  lossy  multiconductor  transmission  line,  con¬ 
sisting  of  current  sources  in  parallel  with  coupled  conductances,  can  be  derived  starting 
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Figure  3.3  Element  stamp  associated  with  the  Thevenin  companion  model  of  Figure  3. 
This  stamp  is  for  the  far  end  of  the  line  only. 
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from  Eqs.  (3. 11a)  and  (3. 1  lb).  These  are  repeated  here  for  convenience  with  the  indepen¬ 
dent  variables  explicidy  shown: 

iCc.tw)  =  (3.30a) 

YQiQi)\ix,co)  =  (3.30b) 

Following  a  development  analogous  to  that  conducted  for  the  modal  voltage  equation,  one 
arrives  at  the  following  equations  that  describe  the  line-end  currents  as  functions  of  the  line- 
end  voltages,  the  reflected  currents  at  the  opposite  end,  and  the  transfer  function  matrices 
Yo(co)  and 

1(0, tu)  =  V(d,G))  +  Yo((u)V(0.ro)  (3.31a) 

l{d,o))  =  -2e"''^“^'V(0,ty)  Yo(£u)V(d,tu)  (3.31b) 

Transforming  these  to  the  time  domain  and  discretizing  them  in  a  manner  similar  to  the 
development  of  the  Thevenin  companion  model  lead  to 

1(0, ffc)  =  -2ri0,tk)  +  h(0,tk)  +  GoV(0,r*)  (3.32a) 

I(d,rJ  =  -2I^(d,rfc)  4-  +  GoV(d,i^)  (3.32b) 


where 


Go  =  lim  Yo(<y) 

<i)— ♦oo 


Yo,(tu)  =  Yo(<u)  -  Go 
/t-i 

r(o,rfc)  =  ^if{d,ti)rid,tk-ti^y)At 
1=0 
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1=0 

k-\ 

I»W.<»)  =  XYo.(<,)V(<i,rrt„,)Ai 

t=0 


The  companion  model  corresponding  to  Eqs.  (3.32a)  and  (3.32b)  is  shown  in  Fig.  3.4, 
and  the  corresponding  eleti -mt  stamp  for  the  far  end  of  the  line  is  given  in  Fig.  3.5.  Again, 
the  topology  of  the  stamp  for  the  near  end  of  the  line  is  the  same  as  that  for  the  far  end. 

3.3.1  Interpretation  of  the  impulse  response  matrices 

It  is  interesting  to  consider  the  physical  interpretation  of  the  impulse  response 
matrices  H(d,r)  and  Zo(r)  encountered  in  this  section.  Referring  to  Fig.  3.6,  a  doubly-infi- 
nite  (n+l)-conductor  transmission  line  with  a  voltage  source  of  value  2, 5(0,  where  5(0  is 
the  impulse  function,  is  shown  inserted  into  line  j  at  x  -  0.  Upon  application  of  the 
impulse  at  r  =  0,  a  voltage  of  value  5(0  is  introduced  on  line  j  at  x  =  0"^,  and  one  of  value 
-5(f)  is  introduced  at  x  =  0~.  The  signal  at  x  =  0"^  will  then  propagate  in  the  positive  x 
direction.  From  the  solution  to  the  multiconductor  transmission  line  equations  given  by 
Eq.  (3.8),  and  making  use  of  the  causality  condition  and  the  inverse  Fourier  transform,  the 
value  of  the  voltage  at  x  =  d  is  easily  found  to  be 

V(d.f)  =  H(d,f)  ♦  [0  •  ■  •  0  5(0  0  ■  •  •  0]  ^  (3.33) 

j 

where  the  impulse  occurs  in  the  yth  row  of  the  vector.  From  Eq.  (3.33)  it  is  een  that 
H,j(d,t),  the  lyth  enay  of  H(d,f),  is  the  voltage  at  x  =  d  of  linr  i  in  the  circuit  of  Fig.  3.6. 

Similarly.,  another  doubly-infinite  f«+l)-conductor  transmission  line  is  shown  in 
Fig.  3  7,  this  time  with  a  current  .source  of  value  25(f)  connected  between  line  j  and  the 
reference  conductor  at  x  =  0.  Upon  application  of  the  impul.se  at  f  =  0,  half  the  current  will 
travel  in  the  positive  x  direction  and  the  other  half  will  flow  in  the  negative  x  direction.  The 
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Figure  3.5  Element  stamp  associated  with  the  Norton  companion  model  of  Figure  3.4. 
This  stamp  is  for  the  far  end  of  the  line  only. 
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Figure  3.6  Physical  interpretation  of  the  impulse  response  Hi^d,t). 
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Figure  3.7  Physical  interpretation  of  the  impulse  response 
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current  flowing  in  the  positive  direction  will  see  an  impedance  equal  to  the  characteristic 
impedance  of  the  line  or  Zo(tu).  Therefore,  in  the  time  domain  the  voltage  at  =  0  is 

V(0,r)  =  Zo(r)  *  [0  •  •  •  0  5(r)  0  •  •  •  0]  ^  (3.34) 

j 

where  again  the  impulse  occurs  in  the  jth  row  of  the  vector.  From  Eq.  (3.34)  it  is  seen  that 
Zoy(r),  the  ijth  entry  of  Zo(t),  is  the  voltage  at  a:  =  0  of  line  i  in  the  circuit  of  Fig.  3.7. 

Since  the  multiconductor  transmission  lines  in  Figs.  3.6  and  3.7  are  infinite  in  each 
direction,  no  reflections  will  occur  at  any  time.  In  a  lossless  system,  the  impulse  response 
is  an  impulse  at  t  =  0  for  each  entry  of  Zo(t)  and  the  sum  of  n  impulses  delayed  in  time  for 
each  entry  of  In  general,  the  delays  are  of  different  length.  Since  there  are  no 

reflections  in  the  responses,  their  durations  are  finite;  if  the  lines  had  been  terminated  by  an 
impedance  other  than  Zq(co),  the  responses  would  exist  for  all  time.  In  the  case  of  lossy 
lines,  each  impulse  response  will  have  a  “tail”  that  is  due  to  the  distortion  of  the  line.  In 
some  sense,  the  length  (in  time)  of  this  tail,  until  the  value  of  the  response  is  negligible,  is  a 
measure  of  the  severity  of  the  distortion  caused  by  the  line.  As  with  the  lossless  case,  the 
impulse  responses  will  still  be  relatively  short  since  no  reflections  occur.  The  shorter 
impulse  responses  make  it  simpler  to  carry  out  the  inverse  Fourier  transform  and  also 
decrease  the  amount  of  past-history  values  that  need  to  be  stored.  This  is  an  inherent 
advantage  of  the  modeling  method  presented  in  this  chapter. 

3.4  Numerical  Considerations 

This  section  contains  a  discussion  of  what  is  involved  numerically  in  calculating 
and  implementing  the  lossy  multiconductor  transmission  line  models  in  a  general  circuit 
analysis  program. 
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T 

The  actual  calculation  of  the  matrix  functions  or  and  Zo(ty)  or  Ygio)) 

is  quite  involved.  The  derivation  from  Eq.  (3.7)  to  Eq.  (3.8)  yields  the  means  of  com¬ 
puting 

^-ncoH  _  E(6))e“^^"^'^E(tu)'^ 

Since  it  need  not  be  considered  further.  Recalling  the  definitions  of 

Zo(<u)  and  Yo(6))  from  Eqs.  (3.10)  and  (3.12): 

Zo(tu)  =  r(<u)-'Z(tu)  Yo(tu)  =  [r(tu)^j-^Y(ty) 

then  making  use  of  Eq.  (3.5)  and  some  matrix  algebra  lead  to 

Zo(tu)  =  E(£0)A(<u)"*E(a))~^Z(C!>) 

Yo(tu)  =  J(n))A(6>)~^J(£t))~^Y(<w) 

However,  because  Yo(ci>)^=  Yo(tu)  and  Y(cu)^=  Y(ty), 

Yo(£0)  =  Y(co)r(£u)“‘  =  Y(G>)E(tu)A(£U)"^E((U)‘^ 

Consequently,  only  the  voltage  eigenvectors  need  to  be  calculated  no  matter  which  trans¬ 
mission  line  model  is  used.  Before  calculating  the  inverse  transform  of  the  variables  above 
by  a  fast  Fourier  transform,  values  must  be  computed  at  discrete  frequencies.  This  is 
complicated  by  the  fact  that,  in  general,  the  eigenvalues  and  eigenvectors  of  Hto)  are  fre¬ 
quency  dependent.  The  numerical  computation  of  the  eigenvalues  and  eigenvectors  of  a 
complex  matrix  is  a  difficult  and  time-consuming  procedure;  however,  because  they  remain 
constant  throughout  a  simulation,  the  impulse  responses  need  only  be  »,alculated  once  for 
each  lossy  multiconductor  transmission  line  of  the  circuit. 

Because  the  transmission  lines  dealt  with  here  are  linear,  the  companion  model  for 
the  lossy  (n-(-l)-conductor  transmission  line  need  only  be  updated  at  each  time  step  of  a 
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simulation  and  not  at  each  Newton-Raphson  iteration.  The  summations  in  Eqs.  (3.26a), 
(3.26b),  (3.29a),  and  (3.29b)  must  be  performed  at  each  time  step.  Over  the  course  of  a 
simulation,  each  summation  constitutes  one-half  of  a  convolution  integral  being  performed 
in  the  time  domain.  These  convolutions  represent  the  most  time-consuming  operation  of 
the  analysis  technique  described  herein.  However,  newer  computers  that  use  vectorized 
arithmetic  units  have  the  capability  of  performing  multiplications  in  parallel.  The  convo¬ 
lution  as  developed  in  this  investigation  is  essentially  a  dot  product  which  is  a  vectorizable 
operation.  Existing  computers  have  the  capability  of  speeding  up  this  operation  by  a  factor 
of  fifty.  In  the  near  future,  a  speedup  of  a  full  two  orders  of  magnitude  may  be  expected. 

The  calculation  of  the  inverse  Fourier  transforms  is  the  most  difficult  numerical 
procedure  involved  in  the  formulation  of  the  time-domain  models.  In  general,  each  entry  of 
ZqCO  or  YqU)  contains  an  impulse  and  a  discontinuity  at  r  =  0.  The  impulse  occurs  because 
Zq(co)  and  Yo(tu)  approach  constant  nonzero  matrices  as  ct)  nears  infinity.  These  impulses 
can  be  extracted  easily  enough  by  shifting  the  frequency  domain  counterpart  of  these  matri¬ 
ces  as  was  done  in  Sec.  3.3.  A  simple  procedure  can  be  used  to  remove  the  discontinuity 
from  the  shifted  response  that  remains.  Since  the  impulse  responses  are  real,  it  follows  that 
Zoi*(<y)  =  Zq^(-co),  for  example,  where  the  asterisk  as  a  superscript  indicates  complex 
conjugate;  then 

2Re{ZQj(£y)}  =  Zq^(q))  +  Zg,  (co)  =  Zq^(co)  +  Zg^I-tu) 

Now  the  inverse  Fourier  transform  of  Zg^I-ft))  is  Zo^(-t);  therefore,  2  RelZo^Itu))  and 
Zos(0  +  Zq^-O  form  a  Fourier  transform  pair.  But  Zo,(0  and  Zg^I-O  do  not  overlap 
because  Zq^IO  is  causal;  as  a  result,  ZqjIO  can  be  recovered  by  taking  the  positive-r  portion 
of  the  Zq/O  +  Zq/-!)  sum.  The  sum  of  these  two  terms  is  real  and  even;  consequently, 
the  discontinuity  at  the  origin  has  been  eliminated.  Similar  reasoning  applies  to  the  com¬ 
putation  of  Yoj(t). 
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Special  problems  arise  whenever  R  or  G  is  zero.  In  many  practical  cases,  a  solu¬ 
tion  of  sufficient  accuracy  may  be  obtained  by  taking  G  to  be  0,  the  zero  matrix.  Upon 
close  examination  of  Eq.  (3.10),  it  is  seen  that  Zo(ci>)  will  then  be  undefined  (infinite) 
whenever  ty  =  0.  However,  Yq(co)  is  the  inverse  of  Zq(co);  consequently,  Yq((0)  will  be 
equal  to  zero  for  co  =  0.  Therefore,  it  would  seem  to  be  advantageous  to  use  the  second 
formulation  of  the  time  domain  model,  the  one  that  results  in  the  Norton  companion  model 
of  Fig.  3.4,  whenever  G  is  zero  or  chosen  to  be  so.  Likewise,  whenever  R  is  zero  at 
to  =  0  and  G  is  not,  it  is  easier  to  use  the  companion  model  of  Fig.  3.2.  If  R  and  G  both 
approach  iiifinity  as  6)  goes  to  zero,  then  in  real  transmission  lines,  both  Yq(co)  and  Zq(co) 
are  also  zero  and  either  method  may  be  used. 

The  inverse  transform  of  the  exponential  propagation  matrix  is  slightly  more 

difficult  to  calculate  than  that  of  the  impedance  or  admittance  matrix.  Imbedded  in  each 
element  of  the  exponential  matrix  are  n  delays  that  correspond  to  the  n  propagation  modes 
of  an  n-line  system.  As  the  frequency  nears  infinity,  the  real  part  of  each  entry  of 
approaches  a  sum  of  n  sinusoids.  As  was  mentioned  in  Sec.  3.3.1,  this  means  that  H(d,t) 
contains  a  sum  of  n  impulse  functions  delayed  in  time  by  the  amount  that  it  takes  the  corre¬ 
sponding  lossless  modes  to  propagate  from  one  end  of  the  line  to  the  other  end.  Once  the 
modal  exponential  propagation  matrix  has  been  recombined  to  the  full-matrix  form, 

there  is  no  simple  procedure  to  extract  these  impulse  functions  as  was  done  for  Yo(ctj)  and 
Zq((o).  In  Yo(<y)  and  Zq(co),  these  impulse  functions  all  occur  at  t  =  0  in  both  the  modal 
and  full-matrix  forms.  The  effect  is  that  goes  to  a  periodic  function  as  co  approaches 

infinity.  Numerically,  the  inverse  fast  Fourier  transform  must  be  of  finite  length;  conse¬ 
quently,  must  be  multiplied  by  some  window  function  before  the  transform  can  be 

calculated.  A  number  of  window  functions  are  available  for  this  purpose.  Among  these 
are  the  rectangular,  triangular,  von  Hann,  and  Hamming  windows.  It  was  found  that  the 
Hamming  window  provided  the  best  results.  The  Hamming  window  w(co)  is  shown  here 


to  be 
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w(cd)  =  0.54  +  0.46  cos  (cDAr)  -nliAt)  <(o<  ;r/(At) 
wico)  =  0  otherwise 

Therefore,  is  calculated  as 

H(d,r)  =  0  f  <  0 

H(d.rt  =  9"‘{2  t  >  0 

where  the  real-and-even  scheme  discussed  earlier  in  the  calculation  of  Zq/O  is  also  used 
above.  The  windowing  results  in  a  minor  spreading  of  the  final  results  that  are  nearly 
imperceptible. 

Each  single  or  multiconductor  transmission  line  represents  an  inherent  delay  in  the 
propagation  of  a  signal  or  signals  during  a  transient  simulation.  This  can  be  <Jeen  in  the 
impulse  response  matrix  H(d,/).  Because  of  this  propagation  delay,  the  value  of  H(d,0)  is 
always  equal  to  the  zero  matrix  for  any  transmission  line  whether  it  is  lossless  or  lossy. 
This  means  that  during  any  transient  simulation,  lumped-element  networks  separated  by 
transmission  lines  are  effectively  decoupled;  that  is,  they  do  not  depend  upon  one  another. 
Consequently,  the  coefficient  matrix  A  can  always  be  reordered  into  a  block  diagonal  form, 
and  each  block  solved  separately  for  each  iteration  of  each  time  step.  If  the  Thevenin 
companion  model  of  Fig.  3.2  is  used,  for  example,  each  half  of  the  equivalent  circuit  can 
be  considered  a  part  of  the  lumped  network  to  which  it  is  connected.  For  large  networks,  a 
significant  savings  of  computer  time  may  result  since  the  number  of  multiplications  re¬ 
quired  to  solve  a  system  of  linear  equations  varies  as  the  cube  of  the  order  of  the  system. 
(If  sparse  matrix  techniques  are  used,  the  number  of  multiplications  for  typical  electrical 
networks  still  varies  at  an  exponential  rate  of  approximately  1.5.)  In  the  typical  circuit 
analysis  program,  the  user  must  number  the  nodes  consecutively  starting  from  one  with  the 
reference  node  being  zero.  In  order  to  ease  the  implementation  of  this  technique  in  the 
program  written  for  this  thesis,  the  user  must  indicate  the  number  of  lumped  element 
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networks  involved  in  the  simulation,  number  the  nodes  in  each  of  these  networks  consecu¬ 
tively,  and  indicate  the  largest  node  value  for  each  lumped-element  network. 

3.5  Conclusions 

In  the  time-domain  analysis  of  nonlinear  circuits,  lossy  multiconductor  transmission 
lines  can  be  represented  at  each  time  step  by  discretized  Thevenin  equivalent  circuits  or 
discretized  Norton  equivalent  circuits.  These  models  contain  only  lumped  elements  so  they 
can  be  easily  implemented  in  a  general  circuit  analysis  program.  The  major  obstacles 
appear  to  be  numerical  in  nature,  but  they  do  not  seem  to  be  severe  enough  to  alter  the 
viability  of  the  technique  presented. 
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CHAPTER  4 

I 

RESULTS  AND  ANALYSIS 


4.1  Introduction 

In  order  to  verify  the  analytical  and  numerical  techniques  developed  in  Chapter  3,  a 
computer  program  compatible  with  the  circuit  analysis  method  of  Chapter  2  was  written 
incorporating  the  lossy  muldconductor  transmission  line  models  [43].  It  is  an  extension  of 
an  earlier  program  written  to  include  lossless  multiconductor  transmission  lines  [44].  The 
new  program  can  be  used  for  the  transient  analysis  of  nonlinear  circuits  containing  any 
number  of  lossy  frequency-dependent  multiconductor  transmission  lines.  It  can  also 
handle  coupled  inductors  and  coupled  capacitors  of  any  order,  as  well  as  junction  diodes 
and  bipolar  transistors.  As  discussed  in  Chapter  2,  the  program  uses  modified  nodal  anal¬ 
ysis  to  formulate  the  circuit  equations,  the  trapezoidal  algorithm  (with  a  fixed  time  step)  for 
discrete  integration,  and  the  Newton-Raphson  algorithm  to  iteratively  solve  systems  of 
nonlinear  equations.  In  this  chapter,  the  program  is  demonstrated  for  several  circuit  con¬ 
figurations. 

4.2  Simulation  Examp’^s 

The  first  example  consists  of  a  single  transmission  line  (with  ground  return)  termi¬ 
nated  by  TI  74S04  (Schottky  TTL)  inverters.  As  shown  in  Fig.  4.1,  one  invener  is  used 
as  the  driver  of  the  line,  and  a  fan-out  of  three  inverters  is  used  to  terminate  the  line.  The 
transmission  line  is  10.0  in.  (0.254  m)  long  with  an  inductance  L  of  307.1  nH/m  and  a 
capacitance  C  of  85.3  pF/m.  The  frequency-dependent  skin  effect  was  modeled  by  a 
resistance  R  =  max {6.8, 0.26^7 /lO^}  Q/in.  =  max {267.7,  10.24^ /lO"^)  f2/m  where/ is 
the  frequency  in  Hertz.  The  conductance  G  of  the  dielectric  was  assumed  to  be  negligible. 
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time  (ns) 


Figure  4.1  Single  lossy  transmission  line  interconnecting  TI 74S04  inverters.  The  output 
of  the  driver  D  is  made  to  switch  from  low  to  high  at  1  ns  and  then  from  high 
to  low  at  22  ns.  For  the  inverters:  rise  time  =  2.2  ns,  fall  time  =  2.2  ns,  Cjn  = 
5.4  pF,  Cout  =  5.0  pF.  For  the  transmission  line:  /?  =  max{6.8, 0.26-^/10^} 
Q/in.  =  max{267.7, 10.24V7/10^}  Q/m,  L  =  307.1  nH/m,  G  =  0.0  S/m,  and 
C  =  85.3  pF/m.  When  losses  are  neglected  the  characteristic  impedance  is  60  Q. 
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These  parameters  are  those  of  a  high-density  interconnect  used  to  transmit  data  in  fast 
computers  [45]. 

For  computational  efficiency,  a  simplified  model  [46]  was  used  to  represent  the 
digital  devices.  When  the  inverter  is  used  as  a  load,  its  input  is  modeled  as  a  parallel 
combination  of  a  constant  capacitance  and  a  voltage-dependent  nonlinear  resistor  as  shown 
in  Fig.  4.2(a).  When  the  inverter  is  used  as  a  driver,  its  output  is  modeled  in  a  manner 
similar  to  that  of  its  input,  except  that  now  the  nonlinear  resistor  is  time  varying.  This 
model  is  shown  in  Fig.  4.2(b).  As  an  example,  the  v-i  curves  of  the  nonlinear  resistors 
used  for  the  TI  74S04  inverter  are  shown  in  Fig.  4.3.  The  input  curve  is  shown  in  Fig. 
4.3^a),  and  the  output  curves  are  shown  in  Fig.  4.3(b).  The  input  curve,  output-low 
cur/e,  and  output-high  curve  are  generated  from  static  measurements  [47].  The  transition 
curves  between  the  two  output  curves  are  generated  by  the  program. 

When  the  driver  is  initially  triggered  by  the  pulse  (see  Fig.  4.1),  the  device  is  oper¬ 
ating  on  the  output-low  curve  of  Fig.  4.3(b).  During  a  period  corresponding  to  the  rise 
time  of  the  inverter,  the  model  passes  through  the  transition  curves  to  the  output-high  cur^e 
where  it  remains  as  long  as  the  input  to  the  driver  is  held  low.  When  the  input  to  the 
inverter  is  switched  back  to  high,  the  output  switches  back  to  the  low  state  via  the  same 
transition  curves  during  the  fall  time  specified  for  the  device.  A  number  of  researchers 
have  had  good  experimental  agreement  using  a  model  of  this  type  [40],[46]. 

The  computer  simulation  of  the  circuit  of  Fig.  4.1  was  first  conducted  using  the 
methods  developed  in  Cliapter  3  to  model  the  transmission  line  and  then  using  the  lumped 
delay  elements  of  Fig.  4.4  to  model  the  line.  Ten  lumped  delay  elements  were  used  to 
approximate  the  line;  the  low-frequency  resistance  was  used  since  frequency  dependence 
cannot  be  modeled  by  the  delay  elements.  A  comparison  of  the  resulting  waveforms  is 
shown  in  Fig.  4.5.  In  Fig.  4.5(a)  the  two  simulations  are  compared  at  the  drive  end,  and 
in  Fig.  4.5(b)  they  are  compared  at  the  load  end.  Except  for  a  small  difference  in  the 
steady-state  values  for  the  output-high  state  of  the  driver,  the  results  are  nearly  identical. 


Inverter  Input  v-i  Curve 
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6.8  n  7.8  nH 


Figure  4.4  Lumped  delay  element:  one  of  ten  cascaded  to  model  the  lossy  transmission 
line  of  Fig.  4.1. 

An  RC  effect  in  the  rise  and  fall  times  that  results  from  the  conductor  losses  can  be  seen  in 
the  figures.  This  is  a  major  problem  facing  users  of  fast  digital  chips. 

The  next  example  consists  of  the  lossless  lines  terminated  by  linear  resistors  and 
capacitors  as  depicted  in  Fig.  4.6(a).  The  inductance  and  capacitance  matrices  and  resulting 
characteristic  impedance  matrix  for  the  coupled  lossless  lines  are  shown  in  Fig.  4.6(b). 
Even  though  the  primary  goal  of  this  thesis  has  been  to  compute  the  transient  responses  of 
lossy  transmission  lines  with  nonlinear  terminations,  this  example  can  serve  as  a  test  for  the 
program  since  the  results  will  be  compared  to  a  Laplace  transform  solution  and  an  experi¬ 
mental  solution  [32].  The  solution  obtained  when  the  method  of  Chapter  3  is  employed  is 
shown  in  Fig.  4.7(a);  the  solution  by  the  method  of  Laplace  transforms  is  shown  in  Fig. 
4.7(b);  and  the  experimental  results  are  shown  in  Fig.  4.7(c).  The  two  numerical  solutions 
are  nearly  identical  as  one  would  expect.  The  experimental  solutions  agree  within  sufficient 
accuracy. 

The  literature  is  sparse  concerning  experiments  that  consist  of  pulse  responses  to 
lossy  multiconductor  transmission  lines.  However,  the  article  by  Gruodis  and  Chang  [30] 
contains  such  an  experiment,  and  the  next  example  compares  simulated  results  using  the 
methods  developed  in  this  investigation  to  their  experimental  results.  The  circuit  used  for 
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Figure  4.5  Comparison  between  simulated  waveforms  for  the  circuit  of  Fig.  4. 1  when  the 
transmission  line  is  modeled  by  the  method  of  Chapter  3  with  a  frequency- 
dependent  resistance  (curves  A)  and  when  it  is  m^eled  by  lumped  delay 
elements  of  the  type  in  Fig.  4.4  (curves  B) 


(a) 
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Figure  4.6  Coupled  lossless  transmission  line  with  linear  terminations  for  numerical  sim¬ 
ulations  and  experimental  measurement:  (a)  circuit  configuration,  (b)  induc¬ 
tance  and  capacitance  matrices  and  resulting  characteristic  impedance  matrix. 


Figure  4,7  Companson  among  two  types  of  simulated  waveforms  and  experimental 
waveforms  for  the  circuit  of  Fig.  4.6:  top,  waveforms  simulated  by  the 
method  of  Chapter  3;  middle,  waveforms  simulated  by  the  Laplace  transform 
method;  bottom,  experimental  waveforms.  For  the  experimental  waveforms, 
the  attenuarion  is  a  factor  of  40  for  the  signal  line  and  10  for  the  sense  line. 
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the  experimental  setup  is  shown  in  Fig.  4.8.  Gruodis  and  Chang  confined  their  experiment 
to  low  frequencies  (around  1  MHz)  in  order  to  minimize  the  effects  of  discontinuities  in¬ 
herent  in  connecting  the  line  to  the  signal  generator  and  terminating  resistors.  The  trans¬ 
mission  line  made  use  of  a  commercial  telecommunication  cable  consisting  of  four  #22- 
AWG  conductors.  One  wire  of  one  of  the  two  twisted  pairs  of  which  the  cable  consisted 
was  used  as  the  signal  line.  The  other  wire  of  the  same  twisted  pair  was  used  as  the  sense 
line.  The  second  twisted  pair  along  with  the  aluminum-coated  mylar  wrapping  and  braided 
shield  formed  the  ground  return.  The  measured  parameters  of  the  transmission  line  are 
shown  in  Table  4.1.  The  experimental  results  are  given  in  Fig.  4.9(a)  and  the  simulated 
results  are  given  in  Fig.  4.9(b).  As  can  easily  be  seen,  all  waveforms  are  very  nearly 
identical  except  for  the  load  end  of  the  signal  line  in  which  the  simulated  results  do  not 
experience  the  slight  rise-time  degradation  that  occurs  in  the  experimental  case.  However, 
the  results  are  well  within  acceptable  limits  of  accuracy. 

Note  that  for  this  simple  example  of  two  symmetrical  lines,  the  eigenvector  matrix 
E  turns  out  to  be  independent  of  frequency: 


consisting  of  only  the  common  and  differential  modes.  The  method  described  herein  is  not 
limited  to  this  case,  and  a  more  general  case  is  considered  in  the  final  example. 

Of  important  concern  in  the  transmission  of  high-speed  pulse  signals  is  the  effect  of 
discontinuities  in  the  multiconductor  lines  used  to  interconnect  digital  circuits.  In  Fig.  4.10 
a  set  of  FAIR  74F04  (Advanced  Schottky  TTL)  inverters  is  shown  interconnected  by  two 
(4+l)-conductor  transmission  lines  joined  by  a  connector.  Figure  4.11  shows  the  micro- 
strip  used  for  the  multiconductor  transmission  lines  along  with  its  inductance  L  and  capac¬ 
itance  C  matrices.  The  L  and  C  matrices  were  calculated  by  the  method  of  Chan  and  Mittra 
[41],  The  microstrip  dielectric  had  a  relative  permittivity  of  4.5.  The  relative  width  of 
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Figure  4-.8  Coupled  lossy  transmission  line  with  linear  terminations.  A  telecommunication 
cable  consisting  of  four  #22-AWG  conductors  in  two  twisted  pairs  is  used  for 
the  transmission  line.  The  signal  line  consists  of  one  wire  of  one  of  the  twisted 
pairs,  and  the  sense  line  consists  of  the  other  wire  of  the  same  twisted  pair. 
The  second  twisted  pair  along  with  the  aluminum-coated  mylar  wrapping  and 
braided  shield  form  the  ground  return.  The  measured  parameters  of  the  trans¬ 
mission  line  are  shown  in  Table  4.1. 
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Table  4. 1  Coupled  transmission  line  parameters  [30]. 


Freq.  (kHz) 

^11 

Resistance  (0/km) 

/?,2  /?22 

0.0 

64.06 

12.83 

64.06 

2.0 

64.01 

12.19 

64.01 

5.0 

63.93 

12.26 

63.93 

10.0 

64.90 

12.45 

64.90 

20.0 

68.29 

13.25 

68.29 

50.0 

77.93 

17.36 

77.93 

100.0 

101.47 

23.59 

101.47 

200.0 

125.72 

26.09 

125.72 

Inductance  (pH/km) 

Freq.  (kHz) 

^12 

0.0 

510.5 

99.5 

510.5 

2.0 

510.5 

99.5 

510.5 

5.0 

498.0 

98.3 

498.0 

10.0 

492.0 

96.4 

492.0 

20.0 

477.7 

94.3 

477.7 

50.0 

451.1 

81.8 

451.1 

100.0 

435.4 

81.8 

435.4 

200.0 

402.2 

68.8 

402.2 

500.0 

360.1 

58.8 

360.1 

1000.0 

340.9 

58.0 

340.9 

Conductance  (S/km) 

Freq.  (kHz) 

G22 

all 

0.0 

0.0 

0.0 

Capacitance  (nF/km) 

Freq.  (kHz) 

Cu 

^12 

^22 

0.0 

86.2 

-14.0 

86.2 

2.0 

74.0 

-12.1 

74.0 

5.0 

74.5 

-12.7 

74.5 

10.0 

75.0 

-14.1 

75.0 

20.0 

75.2 

-15.9 

75.2 

50.0 

75.6 

-14.0 

75.6 

100.0 

75.4 

-15.3 

75.4 

200.0 

76.8 

-14.7 

76.8 

500.0 

76.8 

-14.2 

76.8 

1000.0 

76.6 

-14.6 

76.6 

time  (microseconds) 


Figure  4.9 


Experimental  [30]  and  simulated  waveforms  for  the  circuit  of  Fig.  4.8.  Curve 
A,  one-half  of  drive  signal;  curve  B,  drive  end  of  signal  line;  curve  C,  lr>ad  end 
of  signal  line;  curve  D,  drive  end  of  sense  line;  curve  E,  load  end  of  sense  line. 
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Figure  4.10  Multiconducior  transmission  lines  joined  by  a  connector  and  terminated  by 
FAIR  74F04  inverters.  The  output  of  the  driver  D  is  made  to  switch  from 
low  to  high  at  1  ns  and  then  from  high  to  low  at  30  ns.  For  the  inverters:  rise 
time  =  2.75  ns,  fall  time  =  1.25  ns,  C-^  =  4.6  pF,  Cgm  =  5.0  pF.  The  geom¬ 
etry  of  the  multiconductor  lines  is  depicted  in  Fig.  4.1 1  along  with  their  L  and 
C  matrices.  When  the  lines  are  lossy,  R  =  50  x  1  Q/m,  where  1  is  the  iden¬ 
tity  matrix,  and  G  =  0  S/m.  The  and  Cp  matrices  for  the  connector  are 
given  in  Fig.  4.12. 
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Figure  4. 1 1  Microstrip  used  for  the  multiconductor  transmission  lines  in  the  circuit  of  Fig. 

4.10:  (a)  cross  section,  (b)  inductance  and  capacitance  matrices,  (c)  charac¬ 
teristic  impedance  matrix  when  the  microstrip  is  lossless. 
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Figure  4.12  Inductance  and  capacitance  matrices  for  the  connector  in  the  circuit  of  Fig.  4.10. 

the  conductors,  spacing  between  the  conductors,  and  the  depth  of  the  dielectric  were  taken 
to  be  one.  When  the  lines  were  lossy,  a  diagonal  resistance  matrix  with  non-zero  values  of 
50  QJm  was  chosen.  The  conductance  of  the  dielectric  was  neglected.  This  is  a  reasonable 
assumption  for  the  frequency  range  (up  to  a  fev/  gigahertz)  of  digital  signal  transmission. 
The  connector  was  modeled  by  a  7t-equivalent  circuit  with  inductance  and  capacitance 
C(.  matrices  as  shown  in  Fig.  4.12.  (The  capacitance  matrix  was  split  between  the  two 
halves  of  the  Jt-equivalent  circuit.) 

Figure  4.13  shows  a  comparison  between  simulated  waveforms  for  the  signal  line 
of  the  circuit  of  Fig.  4.10  when  the  transmission  lines  are  lossless  and  when  they  are  lossy; 
Fig.  4.14  shows  the  same  for  the  sense  line  indicated  in  Fig.  4.10.  In  general,  the  losses 
caused  significant  attenuation  and  smoothing  in  both  the  signal  and  sense  lines,  but  only  a 
minor  reduction  of  rise  and  fall  times.  Although  not  obvious  from  these  plots,  the  effect  of 
the  connector  is  the  small  extra  delay  in  the  propagation  time  of  the  signal.  This  should  not 
be  surprising  since  the  model  for  it  is  a  lumped  delay  clement,  albeit  one  that  includes  cou¬ 
pled  elements. 

4.3  Conclusions 


A  number  of  simulation  examples  were  presented  to  exercise  the  lossy  multicon¬ 
ductor  transmission  line  model  in  circuits  that  contain  linear  and  nonlinear  elements.  The 


Figure  4. 13  Comparison  between  simulated  waveforms  for  the  signal  line  of  the  circuit  of 
Fig.  4.10  when  the  transmission  lines  are  lossless  and  when  they  are  lossy. 


Figure  4.14  Comparison  between  simulated  waveforms  for  the  sense  line  of  the  circuit  of 
Fig.  4.10  when  the  transmission  lines  are  lossless  and  when  they  are  lossy. 
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examples  highlighted  the  interconnection  of  high-speed  digital  circuits  by  lossy  multicon¬ 
ductor  transmission  lines  with  and  without  discontinuities.  Simulation  results  showed 
good  agreement  with  experimental  results  for  both  lossless  and  lossy  multiconductor 
transmission  lines. 
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CHAPTER  5 

CONCLUSIONS 

Using  the  frequency-dependent  transmission  line  parameters,  two  time-domain 
models  were  developed  for  lossy  multiconductor  transmissions  lines.  It  was  shown  that 
the  endpoints  of  a  lossy  multiconductor  line  could  be  represented  at  each  time  step  by  dis¬ 
cretized  Thevenin  or  Norton  equivalent  circuits.  Because  these  models  contain  only 
lumped  elements,  they  can  be  easily  implemented  in  a  general  circuit  analysis  program  for 
simulating  the  transient  responses  of  nonlinear  circuits.  A  number  of  transient  simulations 
were  conducted  in  order  to  test  the  viability  of  the  program.  The  program  was  found  to 
give  good  results  and  to  be  of  a  very  general  nature. 

Of  particular  interest  was  the  interconnection  of  high-speed  digital  circuits  by  lossy 
multiconductor  transmission  lines  containing  discontinuities.  The  technique  was  demon¬ 
strated  for  the  most  general  case  of  lossy  multiconductor  transmission  lines  with  nonlinear 
terminations  in  which  the  modal  transformation  matrices  were  necessarily  frequency 
dependent.  This  case  has  been  avoided  very  adroitly  in  the  literature.  Results  have  only 
been  shown  for  the  simplest  case  of  two  symmetric  conductors  over  a  ground  plane  [30]- 
[32].  In  this  particular  case,  the  modal  transformation  matrices  result  in  an  even  and  odd 
mode  no  matter  what  the  frequency;  therefore,  the  most  general  case  is  avoided.  Also,  the 
capability  of  incorporating  discontinuities  into  the  transmission  line  structure  was  presented 
for  the  first  time. 

A  unique  feature  of  the  analysis  procedure  developed  in  this  study  was  the  exclu¬ 
sive  use  of  infinite-line  impulse  responses  in  the  formulation  of  the  time-domain  models. 
Because  infinite  lines  are  matched,  there  were  no  reflections  in  the  impulse  responses.  The 
result  is  that  these  impulse  responses  are  necessarily  shorter  than  those  for  the  quasi- 
matched  case,  such  as  the  modified  short-circuit  parameters  of  Djordjevic  et  al.  [31]  and  the 
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scattering  parameters  of  Schutt-Aine  [32].  Shorter  impulse  responses  reduce  the  problem 
of  aliasing  encountered  when  taking  the  inverse  transform  and  decrease  the  amount  of 
computation  that  needs  to  be  performed  in  the  time  domain.  Techniques  were  presented  to 
facilitate  the  implementation  of  the  fast  Fourier  transform  which  is  utilized  to  convert 
frequency-domain  transfer  functions  to  time-domain  impulse  responses. 

The  most  serious  drawback  of  the  proposed  technique  is  that  direct  integration  of 
the  convolution  integrals  is  carried  out  in  the  time  domain.  Since  the  number  of  multiplica- 
tions  that  must  be  performed  is  of  the  order  of  N  ,  where  N  is  the  number  of  steps  in  the 
simulation,  this  is  indeed  a  very  serious  consideration.  However,  newer  computers  that 
use  vectorized  arithmetic  units  have  the  capability  of  performing  multiplications  in  parallel. 
The  convolution  as  developed  in  this  investigation  is  essentially  a  dot  product  which  is  a 
vectorizable  operation.  Existing  computers  have  the  capability  of  speeding  up  this  opera¬ 
tion  by  a  factor  of  fifty.  In  the  near  future,  a  speedup  of  a  full  two  orders  of  magnitude  can 
be  expected. 

There  appears  to  be  two  possible  solutions  to  the  problem  mentioned  above.  In  the 
first  of  these,  the  impulse  responses  are  approximated  by  sums  of  exponentials;  then,  the 
convolutions  can  be  accomplished  recursively.  In  the  second,  the  transfer  functions  are 
approximated  by  rational  polynomials  directly  in  the  frequency  domain  so  that  state-variable 
techniques  may  be  used.  In  both  of  these  methods,  the  amount  of  computation  needed  to 
calculate  the  convolutions  is  greatly  reduced  and  the  need  to  store  past-history  values  of 
voltage  and  current  eliminated.  The  second  procedure  has  been  attempted  by  Gruodis  and 
Chang  [30]  in  a  method  similar  to  the  one  presented  herein.  They  use  the  diagonal  form  of 
the  exponential  propagation  matrix.  The  delays  are  extracted  from  the  exponential  func¬ 
tions,  and  modal  signals  are  propagated  along  the  transmission  line  by  shifting  them  at  each 
time  step  as  is  standard  for  lossless-line  simulations.  The  waves  are  then  shaped  at  the 
opposite  end  of  the  line  by  the  remaining  ponions  of  the  propagation  functions.  However, 
in  converting  between  modal  variables  and  line  variables,  they  neglect  the  fact  that  the 
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transformation  matrices  are  frequency-dependent.  A  physical  interpretation  of  these  trans¬ 
formation  functions  is  not  possible;  consequently,  their  approximation  is  necessarily  more 
difficult  since  they  may  be  noncausal. 

The  inclusion  of  the  lossy  multiconductor  transmission  line  models  into  an  existing 
circuit  simulator  such  as  SPICE  is  the  most  obvious  future  need.  Current  computer-aided 
design  tools  for  nonlinear  circuits  have  the  ability  to  model  single  lossless  lines  or  use  RC 
ladders  to  model  single  resistive  lines.  They  cannot  accurately  portray  the  multiline  cou¬ 
pling  and  frequency-dependent  losses  that  arise  when  high-speed  signals  propagate 
between  digital  circuits.  The  method  presented  in  this  study  is  currently  limited  to  uniform 
transmission  lines.  Tapered  lines  are  also  used  by  the  computer  and  telecommunication 
industries;  as  a  result,  the  inclusion  of  nonuniform  lines  is  also  of  importance.  In  addition, 
the  inclusion  of  nonlinear  lines  may  be  needed  eventually,  but  a  totally  new  procedure 
would  have  to  be  devised  since  the  Fourier  transform  could  no  longer  be  used.  This  will 
probably  not  be  necessary  for  digital  circuit  interconnections  in  the  near  future. 
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/lA.t/Kiirr  — This  work  presents  a  new  approach  for  the  timeMlomain 
simulation  of  transients  on  a  ilispersire  and  lossy  transmission  line 
terminated  with  active  devices.  The  method  combines  the  scattering  matrix 
of  an  arbitrary  line  and  the  nonlinear  causal  impedance  functions  at  the 
loads  to  derive  expressions  for  the  signals  at  the  near  and  far  ends. 

The  problems  of  line  losses,  dispersion,  and  nonlinearities  are  first 
investigated.  A  time-domain  formulation  is  then  proposed  using  the  scatter¬ 
ing  matrix  representation.  The  algorithm  assumes  that  dispersion  and  loss 
models  for  the  transmission  lines  are  available  and  that  the  frequency 
dependence  is  known.  Large-signal  equivalent  circuits  for  the  terminations 
are  assumed  to  be  given.  Experimental  and  computer-simulated  results  are 
compared  for  the  lossless  dispersionless  case,  and  the  effects  of  tosses  and 
dispersion  are  predicted. 

I.  Introduction 

N  TODAY’S  MANY  applications  of  integrated  circuits 
and  printed  circuit  boards,  transmission  lines  and  inter¬ 
connections  play  an  instrumental  role  at  virtually  every 
level  of  integration.  With  the  design  of  fast  devices  having 
switching  times  in  the  picosecond  range,  transmitting  data 
at  high  megabaud  rates  has  become  very  commonplace  in 
modem  digital  computers  and  switching  networks  used  for 
telecommunication.  Signal  delays  and  rise  times  are  more 
and  more  limited  by  interconnection  lengths  rather  than 
by  device  speed  and  represent  a  potential  obstacle  to  the 
ultimate  scaling  on  VLSI  technology.  In  recent  years, 
modeling  interconnections  has  become  a  major  focus  of 
interest  in  the  implementation  of  digital  and  microwave 
circuits.  Shorter  rise  and  fall  times  as  well  as  higher 
frequency  signals  have  compelled  most  transmission  lines 
to  operate  within  ranges  where  dispersion  is  no  longer 
negligible.  Skin  effect  and  losses  contribute  to  signal  cor¬ 
ruption  leading  to  waveform  attenuation  as  well  as  pulse 
rise  and  fall  time  degradations.  In  wafer-scale  integration, 
these  losses  can  become  very  significant  and  may  lead  to 
an  RC  type  behavior  of  the  lines.  Finally,  in  the  case  of 
multiconductor  lines,  cross-coupling  between  neighboring 
lines  may  increase  the  level  of  distortion  in  excited  lines 
which  can  initiate  false  signals  in  nonexcited  lines. 

The  implementation  of  a  high-density-compatible  pack¬ 
aging  scheme  is  essential  for  the  design  of  high-speed 
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digital  systems  such  as  gallium  arsenide  integrated  circuits. 
For  microwave  or  digital  applications,  printed  circuit 
boards,  chip  carriers,  and  modeling  of  these  networks 
represent  the  first  step  toward  implementing  reliable  de¬ 
sign  guidelines.  A  complete  CAD  tool  for  studying  these 
effects  would  require  a  frequency-domain  characterization 
of  the  transmission  line  with  higher-order  modes  included 
to  account  for  dispersion.  Numerous  authors  have  investi¬ 
gated  the  properties  of  microstrip  lines  at  high  frequencies 
and  derived  expressions  relating  the  propagation  character¬ 
istics  to  frequency  (11-[10].  Full-wave  and  simplified  mod¬ 
els  have  been  proposed  to  describe  these  effects  and  to 
derive  the  frequency  dependence  of  the  characteristic  im¬ 
pedance  and  the  propagation  constant.  Other  geometries, 
such  as  stripline,  buried  microstrip,  and  coplanar,  have 
thus  far  received  less  attention  but  obey  the  same  restric¬ 
tions  imposed  on  the  electrical  performance  of  microstrip 
at  microwave  frequencies. 

Of  equal  importance  is  the  analysis  of  a  high-speed  or 
high-frequency  signal  propagating  on  a  dispersive  and 
lossy  transmission  line.  Such  an  analysis  requires  a  com¬ 
plete  and  accurate  frequency  characterization  of  the  struc¬ 
ture  of  interest  and,  for  practicality,  must  implement  the 
nonlinear  and  time-changing  behaviors  of  the  termina¬ 
tions,  which  are  transistors,  logic  gates,  or  other  types  of 
active  devices.  Several  investigators  have  attempted  to  set 
up  analytical  models  describing  wave  propagation  in  such 
systems.  Solutions  for  lossless  lines  with  arbitrary  termina¬ 
tions  were  obtained  by  Mohammadian  el  al.  (12]  using  a 
forward  and  backward  wave  approach.  Veghte  and  Balanis 

[13]  have  analyzed  the  distortion  of  a  pulse  due  to  dis¬ 
persion  along  a  microstrip  transmission  line.  Caniggia 

(14)  combined  macromodels  for  transmission  lines  and 
terminations,  and  Djordjevic  et  al.  [11]  used  a  Green's 
function  approach  to  simulate  the  time-domain  transient 
on  a  multiconductor  array  with  nonlinear  terminations. 

In  this  study,  a  combined  frequency-domain,  time- 
domain  approach  is  used  to  formulate  the  propagation 
equations  on  a  dispersive  and  lossy  line  with  nonlinear  be¬ 
havior  at  the  terminations.  The  novelty  of  the  method 
resides  in  the  formulation,  which  separates  the  linear  ex¬ 
pressions  for  the  transmission  line  from  the  nonlinear 
expressions  for  the  terminations  by  means  of  scattering 
parameters  and  the  use  of  an  auxiliary  reference  imped- 
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ance.  A  time-domain  flow-graph  representation  of  the 
solution  is  also  derived.  In  the  lossless  case,  the  solution 
reduces  to  very  simple  expressions  which  greatly  increase 
computational  efficiency. 

II.  Formulation 

Consider  an  arbitrary  transmission  line  with  arbitrary 
loads  at  both  ends  (sec  Fig.  1).  The  differential  equations 
relating  the  voltage  V  and  the  current  I  along  the  line  are 
expressed  by 


dV 

01 

dx 

(la) 

di 

ov 

(lb) 

where  Lq,  Cq,  Rq,  and  Gq  are  the  inductance,  capacitance, 
resistance,  and  conductance  per  unit  length,  respectively. 
The  solutions  for  time-harmonic  excitation  are  usually 
written  in  the  frequency  domain  as  (u  =  2?r/  is  angular 
frequency) 

Viu,x)  ^  Ae-''^  +  Be*^'  (2a) 

1 

/(u.x)  Se*’’-']  (2b) 

••vhere 

; -  I  Rn-^  yw'f-o 

V  Gq  +  juCo 

(3) 

Here  y  .ind  Zg  are  complex,  leading  to  an  attenuation  of 
the  signal  as  it  propagates  through  the  medium.  If  the 
terminations  are  linear  and  time-invariant  (i.e.,  Z,(t)  and 
Z-(/)  are  constant  with  time),  the  coefficients  and  B 
can  be  determined  by  matching  boundary  conditions  at 
.r  =  0  and  x  =  /;  next,  an  inverse  Fourier  transform  ap¬ 
proach  can  be  used  to  solve  for  the  time-domain  solution. 
On  the  other  hand,  if  the  terminations  are  nonlinear  or 
time-changing,  then  the  boundary  conditions  must  be  for¬ 
mulated  in  the  time  domain  as 

F,(r)-K(r.0)+Z,(/)/(r.0)  (4) 

Z:(r)/(f./)  =  K(r,/)  (5) 

where  y,(r)  is  the  source  voltage.  Zjfr)  and  Z,(r)  indicate 
the  time  variations  of  the  source  and  load  terminations, 
respectively.  For  any  time  greater  than  r,  Z,(r)  and  Zjfr) 
are  not  known,  since  they  depend  on  the  voltage  .and 
current  solutions  at  time  t.  In  fact,  the  evaluation  of  Z,(i) 
and  Zgfr)  may  involve  several  iterations  that  involve  solv¬ 
ing  the  terminal  network  equations  with  trial  values  until 
convergence  to  the  true  impedance  values.  Transforming 
conditions  (4)  and  (5)  Into  the  frequency  domain  is  inap¬ 
propriate,  and  a  time-domain  formulation  is  thus  neces¬ 
sary.  Likewise,  (2)  and  (3)  cannot  conveniently  be  analyti¬ 
cally  inverted  into  the  time  domain  and  constrained  to 
satisfy  (4)  and  (5).  This  limitation  arises  not  only  because 
of  the  dispersive  and  frequency-dependent  characteristics 
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Fig.  1.  Transmission  line  with  nonlinear  terminations  and  source 

generator. 


of  the  line,  but  also  because  the  evaluation  of  the  coeffi¬ 
cients  A  and  B  requires  an  a  priori  knowledge  of  the 
time-dependent  load  functions.  A  formulation  in  which  the 
causality  of  the  boundary  conditions  is  implied  thus  be¬ 
comes  necessary.  The  use  of  scattering  parameters  allows 
one  to  define  the  properties  of  the  transmission  line  inde¬ 
pendently  from  those  of  the  terminations;  consequently, 
by  properly  combining  load  and  line  relations,  a  simple 
expression  for  the  solutions  can  be  derived. 


HI.  Scattering  Parameter  Formulation 


Any  linear  two-port  network  can  be  described  as  a  set  of 
scattering  parameters  (S  parameters)  which  relate  incident 
and  reflected  voltage  waves.  These  waves  are  variables 
which  depend  on  the  total  voltages  and  currents  at 
the  two-port.  If  ideal  (lossless  and  dispersionless)  transmis¬ 
sion  lines  of  known  characteristic  (reference)  impedance 
are  connected  to  both  ports  of  a  linear  network,  then 
the  voltage  waves  on  the  reference  lines  (see  Fig.  2) 
flj,  h;,  a,,  hi  are  defined  as  the  incident  and  reflected 
waves  from  port  1  and  port  2,  respectively.  The  scattering 
parameters  are  then  known  to  satisfy  the  frequency-domain 
relation 


=  (6) 

h,  ”  5i,fl, -t- 5^,0,.  (7) 

S,,  and  5ii  are_  regarded  as  scattering  reflection  coeffi¬ 
cients,  whereas  S,,  and  are  the  scattering  transmission 
coefficients  of  the  network.  The  total  voltages  in  ports  1 
and  2  are,  respectively,  given  by 

F,  =  a,-l-h,  (8) 

F,-a,  +  h,  (9) 

and  the  expressions  for  the  currents  are 


/,= 


'ref 


Zref 


(10) 


(11) 


A  relation  between  the  propagation  characteristics  of  a 
transmission  line  and  the  associated  scattering  parameters 
can  then  be  easily  derived.  For  a  single  mode  of  propa¬ 
gation,  it  can  be  shown  that  for  a  given  transmission  line 
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(1-P~)« 
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Fig.  2.  Transmission  line  and  frequency- domain  flow-graph  representa¬ 
tion  using  scattering  parameters. 
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Fig.  3.  Time-domain  circuit  and  flow-graph  representation  of  arbitrary 
transm'ssion  line  (TL)  with  nonlinear  terminations  at  time  t.  The  sign  * 
indicates  a  convolution  between  time-domain  scattering  parameters 
and  the  independent  voltage  waves  a^  and  a-  as  per  (14)  and  (15). 


where  a(u)  and  Zq(o))  account  for  the  dispersive  and 
lossy  behaviors  of  the  line.  Zq(u)  is  the  characteristic 
impedance  of  the  line  to  be  analj'zed  and  must  be  dis¬ 
tinguished  from  Z,,f,  the  characteristic  impedance  of  the 
reference  lines  which  support  the  voltage  waves  associated 
with  the  5  parameters.  The  scattering  parameters  of  a 
transmission  line  depend  only  on  its  electrical  characteris¬ 
tics  and  are  not  influenced  by  the  source  and  load  voltages 
at  the  terminations;  however,  the  overall  response  of  the 
system  is  a  combination  of  line  and  termination  responses 
and  can  be  obtained  by  cascading  the  various  sections  of 
the  network.  We  can  then  write  the  time-domain  equations 
relating  the  voltage  waves  of  an  arbitrary  line  terminated 
with  nonlinear  loads  (Fig.  3).  We  get  (subscripts  I  and  2 
refer  to  near  end  and  far  end,  respectively) 

f>,(0=5H(f)*a,(f)  +  5,,(f)*o,(f)  (14) 

ft:(0=52,(f)*ai(f)  +  S,,(f)*C:(f)  (15) 

where  *  indicates  a  convolution  in  the  time  domain.  The 
scattering  parameters  5,i(f),  S,,(f),  S2,(/),  Szzff )  are  the 
inverse  transforms  of  the  frequency-domain  S  parameters 
ar.d  can  be  viewed  as  Green’s  functions  associated  with  the 
time-domain  response  of  the  transmission  line,  due  to  a 
single  frequency  source  at  the  terminations.  The  load 
conditions  at  the  near  and  far  ends  are  now  directly 
expressed  in  the  time  domain  by  looking  at  the  flow-graph 
representation  of  the  system  (see  Fig.  3) 

«i(0  =  r,(f)/>,(f)-4-r,(/)i;,(f)  (16) 

O,(f)  =  r,(f)/i2(0  +  r,(f)g,(f)  (17) 

in  which  r,(f ).  r,(f ),  r,(/),  Tjff)  are  the  reflection  and 
transmission  coefficients  associated  with  near  and  far  ends. 


Z,(0  +  2«f 


r,(/) 


Z|(f)+2„f 


(18) 


In  (14)  and  (15),  each  of  the  convolution  terms  can  be 
expressed  as 

~  (19) 

Since  the  algorithm  to  be  derived  must  be  amenable  to 
computer  usage,  it  is  desirable  to  discretize  (19)  and  isolate 
aj(/)  in  a  manner  analogous  to  [11,  eqs.  (12)  and  (13)]; 

•5//(0*fly(f)  “  E  ■J//f-T)ay(T)AT 

T-l 

•5/>(0*a>(0  “•9,y(0)fly(OAT-f-  Z 

T-l 

(20) 


or 

5.v(/)*«,(0  “■?,'(0)o,(/)+  //,,(0  (21) 

where  At  is  the  time  step  and  Sf  (0)  =•  S,y(0)AT.  //,y(t)  = 
T.rZ.\Sij(i  -  r)aj(r)^r  represents  the  history  of  the  line 
and  depends  on  information  up  to  time  t-l.  Causality 
insures  that  the  Oj's  are  known  tor  r  <f,  which  allows  the 
use  of  this  information  for  the  determination  of  the  a^’s  at 
T  =  /.  We  first  substitute  (21)  into  (14)  and  (15)  and  obtain 

6,(r)-S,',(0)ai(/)  +  5;2(0)o2(r)+f/„(r)  + //„(;) 

(22) 

62(f)  =5':f,(0)o,(/)  +  .yA(0)o2(/)  +  //2,(/)  +  //22(/). 

(23) 


Combining  the  above  equations  with  those  for  the  forward 
waves  (16)  and  (17),  one  gets 

^  ^  [l-r2(r)5f2(0)][7',(r)g,(f)-f  r,(r)A/,(r)] 


r,(r)s;2(o)[r2(t)g2(0  +  r2(t)A/2(t)] 

^ - M? -  <-■'> 

,  ^  [i-r,(0^;,(o)l(n(,)sj(0+r,(,)M,(01 

- MO - 

r2(f)A’:f.(0)[r,(r)g,(O  +  r.(OA/,(/)] 

A(r) 

A(o-[i-r,(t)s;,(o)][E-r,52'2(o)] 

-r,(r)S,'2(0)r2(t)52',(0)  (26) 


where  ,!/,(/)  =  //„(() -h //^(t),  and  A/,(t )  =  //,,(/) + 
HtjU).  The  variables  h,(/)  and  5,(0  ^re  recovered  using 
(22)  and  (23),  and  the  total  voltages  at  ports  1  and  2.  by 
using  (8)  and  (9).  Fig.  4  shows  the  flow-graph  representa¬ 
tion  for  the  transmission  line  at  time  /  in  which  the 
memory  of  the  complete  network  has  been  included  in  the 
terms  A/,(t)  and  A/,(/).  The  independent  terms  are  g,(/) 
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Fig.  4.  Equivalent  flow-graph  representation  of  arbitrary  transmission 
line  with  nonlinear  lime-varying  terminations  at  time  t.  The  history  of 
the  networic  is  included  via  the  memory  variables  ^/,(r)  and  Mtfl). 
Note  that  this  representation  involves  only  multiplications  in  the  time 
domaia 


and  ^2(0-  Wi(0  ni'd  fire  also  independent  and 

contain  the  information  pertinent  to  the  history  of  the  line. 
Numerical  efficiency  is  of  practical  importance  for  the 
simulation  since  it  determines  the  speed  of  the  computa¬ 
tions  involved.  Expressions  (24)-(26)  can  be  further  re¬ 
duced  by  observing  that  for  transmission  lines  with  non¬ 
zero  length,  S{^(0)  and  S^iiO)  must  vanish,  since  a  finite 
duration  is  required  for  an  arbitrary  signal  to  propagate 
through  the  line.  The  above  relations  then  become 

,  .  [l-r,(/)5i(0)|[7-,(.)j,(r)+r,(,)M,(.)i 

- m7) - ^ - 

(27a) 


,  ,  [i-rt(f)5(i(o)][r2(/)g2(/)+r2(/)M2(0] 

- - MO - 

(27b) 

&,(/)-S(,(0)u,(f)  +  A/i(f)  (27c) 

^i(0  =M2(0)«:(0  +  A/2(0  (27d) 


with 

ti(r)  =  [l-rj(f)5(,(0)l[l-r2(/)M2(0)l-  (28) 

Computational  limitations  in  these  expressions  are  es¬ 
sentially  determined  by  M^{t)  and  M^O),  which  con¬ 
tain  the  history  of  the  network  and  involve  (he  voltage 
wave  solutions  from  previous  time  steps.  Tn  the  case  where 
losses  and  dispersion  are  neglected,  the  frequency-domain 
scattering  parameters  associated  with  the  transmission  line 
become 

(l-e--'''“'/'’)p 

p2g-2yw/A'  (^9) 

5,2(01)  (_-p2^-2;WA  (30) 

where 


^0  ~ 

^0  ^ref 


(31) 


Since  Zg,  the  characteristic  impedance  of  the  line,  is  con¬ 
stant  with  time,  and  since  Z„f,  the  reference  impedance,  is 
arbitrary,  one  can  choose  Z„f  *  Zq,  which  leads  to  p  =  0 
and 


Sn(o>)- 522(01) -0  (32) 

(33) 


Therefore,  the  time-domain  Green's  functions  associated 


with  the  scattering  parameters  are 

•5n(0 -*522(0-0  (34) 

5,2(0 -.52i(0-8('-^)  (35) 

M,(t)~a,(t-I/v)  (36) 

W2(/)-ni(f-//0  (37) 

MO-1.  (38) 


We  then  obtain  simple  expressions  for  the  forward  and 


backward  waves: 

"i(0  -7’i(/)g,(/)  +  r,(/)a2(^~/A) 

(39) 

0 20)~ '^10)82(1)+ ^2(1) 

(40) 

*i(0-ai(f-V") 

(41) 

^>2(0 -ai(f-^A). 

(42) 

The  advantage  of  the  above  expressions  lies  in  their 
computational  efficiency,  since  only  a  search  is  involved  in 
the  evaluation  of  the  history  of  the  network,  and  no 

summations  of  previously  calculated  terms  are  needed. 

IV.  Models  for  Terminations  and  Devices 

Thus  far,  this  study  concentrated  on  simulating  the 
time-domain  transient  response  for  arbitrary  transmission 
lines  terminated  with  nonlinear  devices.  In  this  section,  wc 
examine  the  nature  of  the  terminations  and  the  manner  in 
which  they  are  to  be  represented  in  a  form  consistent  with 
the  relations  derivecl.  Several  techniques  are  available  that 
convert  reactive  elements  and  nonlinear  devices  to  time- 
varying  causal  resistances  as  well  as  voltage  or  current 
sources.  We  briefly  overview  two  of  these  techniques, 
namely,  the  trapezoidal  algorithm  and  the  Newton- 
Raphson  (NR)  algorithm.  A  more  detailed  development 
can  be  found  in  [15].  Since  formulation  and  solution  are  in 
the  time  domain,  every  element  must  have  an  equivalent 
network  in  the  time  domain.  The  trapezoidal  algorithm  is  a 
numerical  integration  algorithm,  and  its  use  in  repre¬ 
senting  reactive  elements  in  the  time  domain  is  illustrated 
below.  Consider  a  capacitor  C,  with  a  current-voltage 
relation  given  by 

dV 

I.C-.  (43) 

If  we  discretize  the  time  variable  by  choosing  a  time  step 
h,  then  the  voltage  at  time  i,+i“(n+l)A  can  be 
approximated  in  terms  of  variables  at  r »  nit  as 

(44) 

where  the  superscript  '  indicates  a  derivative  with  respect 
to  time.  Making  use  of  (43),  we  get 


Equation  (45)  can  then  be  represented  by  the  equivalent 


schutt-aine  and  mittea:  scatterino  parameter  transient  analysis 


S33 


Fig.  5.  Linear  capacitor  and  equivalent  lime-domain  irapezoidal  al¬ 
gorithm  representation. 


Fig.  6.  Geometrical  and  circuit  interpretation  of  the  Newton- Raphson 
algorithm  at  the  k  th  iteration. 


linear  one-port  model  at  time  f„+i=’(rt  +  l)A  (see  Fig.  5) 
2C  1C 

with  conductance  G  =  —  and  current  source  J„  =>  —  F, 

-f  /„.  An  analogous  derivation  can  be  performed  for  induc¬ 
tances  as  well.  The  efficiency  of  the  method  depends  on  a 
proper  choice  of  the  time  step  which  determines  the  stabil¬ 
ity  of  the  numerical  solution.  Several  other  techniques  of 
numerical  integration  are  also  available  which  offer  greater 
stability  at  the  price  of  numerical  complexity  [15]. 

Nonlinear  elements  such  as  diodes  and  transistors  must 
be  reduced  to  equivalent  networks  with  linear  elements  at 
time  I.  If  the  nonlinear  current  voltage  relation  for  these 
elements  is  known,  then  an  iterative  scheme  such  as  the 
Newton-Raphson  algorithm  can  be  used  to  seek  a  solu¬ 
tion.  The  circuit  representation  of  the  Newton-Raphson 
technique  is  illustrated  in  Fig.  6.  At  a  particular  time,  a 
guess  value  for  the  voltage  is  chosen  to  which  a  current  is 
associated  via  the  /-k'  relations  that  determine  the  operat¬ 
ing  point  Q.  The  next  guess  is  then  related  to  the  previous 
one  by 


(46) 


At  each  iteration  step,  the  resulting  equivalent  circuit 
is  composed  of  a  linear  conductance  of  value  G,.  =  dl/dV 
at  and  a  current  source  with  value  given  by  » 
G^^^.  Solving  the  combined  transmission  line 
Newton-Raphson  equivalent  circuit  problem  at  each  itera¬ 
tion  step  will  lead  to  the  actual  representation  of  the 
termination  at  time  t. 

Once  linearization  and  discretization  are  performed, 
time-domain  values  are  available  for  the  equivalent  resis¬ 
tances  or  generators.  These  expressions  are  causal,  since 
their  values  at  any  time  t  depend  on  the  history  of  the 
network  which  renders  impossible  an  a  priori  knowledge 
of  the  time  variations  of  ihe  termination  impedance.  Non¬ 
linear  complex  elements  can  be  handled  by  first  using  the 


NR  scheme  for  linearization  at  a  given  time,  then  stepping 
in  the  time  domain  while  replacing  linear  complex  ele¬ 
ments  by  time-varying  resistances  arid  generators. 

V.  Appucations  and  Practical  Considerations 

Many  applications  in  microwave  and  digital  communi¬ 
cations  require  the  use  of  transmission  lines  terminated 
with  nonlinear  devices.  Distortion  and  noise  arise  when  the 
terminations  are  not  matched  to  the  line  impedance.  More¬ 
over,  if  losses  and  dispersion  are  present  in  the  line, 
attenuation  and  time  delay  come  into  account.  The  combi¬ 
nation  of  these  effects  needs  to  be  modeled  and  simulated 
on  a  reliable  computer-aided  design  (CAD)  tool.  Simula¬ 
tions  of  pulse  propagation  through  lossy  transmission  lines 
terminated  with  active  devices  can  be  very  useful  in  pre¬ 
dicting  signal  distortion,  attenuation,  rise  and  fall  time 
degradation  which  occur  along  the  transmission  path.  The 
necessary  information  for  such  a  tool  are  the  line  parame¬ 
ters  over  a  wide  frequency  range  and  complete  device 
characteristics.  The  line  parameters  can  be  found  from  a 
frequency-domain  full-wave  dispersion  analysis  which  in¬ 
cludes  the  effects  of  losses.  Device  data  are  usually  ob¬ 
tained  from  the  current-voltage  characteristics  provided 
by  the  manufacturers. 

Several  computer  simulation  programs  were  developed 
to  simulate  waveform  distortion  on  various  combinations 
of  transmission  lines  and  terminations.  First,  a  lossless 
stripline  structure  terminated  with  advanced  Schottky  (AS) 
TTL  inverters  was  studied  (see  Fig.  7).  The  length  of  the 
line  was  50  in  (1.27  m)  with  Zo  =  73  Q,  d  =  0.142  m/ns. 
The  line  was  excited  by  the  output  of  an  AS04  inverter 
(driver).  The  receiving  end  of  the  line  was  connected  to  the 
inputs  of  eight  AS24C  inverters.  For  computational  ef¬ 
ficiency,  a  simple  dynamic  nonlinear  equivalent  circuit  was 
used  to  model  the  output  of  the  driver.  The  network 
consisted  of  a  voltage  generator  in  series  with  a  voltage- 
dependent  resistor.  An  8  pF  capacitor  was  placed  in  paral¬ 
lel  with  the  combination  to  model  rise  and  fall  time 
degradation  and  RC  time  delays  (see  Fig.  7(b)).  The  gener¬ 
ator  provided  a  pulse  with  a  magnitude  of  4.2  V.  The 
voltage  dependence  of  the  resistance  is  shown  in  Fig.  7(b). 
The  quick  jump  in  resistance  is  used  to  model  the  cutoff 
point  of  one  output  transistor  of  the  TTL  inverter.  The 
input  of  each  AS240  inverter  was  modeled  as  a  reverse-bi¬ 
ased  Schottky  diode  in  parallel  with  an  8  pF  capacitor  (see 
Fig.  7(c)).  The  simulation  process  involved  first  choosing  a 
time  interval  and  then  stepping  in  time  and  determining 
the  voltage  variables  as  per  (39)-(42).  The  trapezoidal 
scheme  was  used  to  convert  capacitors  to  linear  sources 
and  resistances.  Likewise,  the  Newton-Raphson  algorithm 
made  it  possible  to  convert  the  Schottky  diodes  into  a 
linearized  equivalent  circuit.  The  time  steps  were  found  to 
have  no  significant  effects  on  the  accuracy  of  the  solutions 
due  to  the  good  stability  properties  of  the  trapezoidal 
algorithm.  Experimental  results  are  compared  with  the 
simulations  in  Fig.  8(a)  and  (b).  Minor  discrepancies  were 
attributed  to  pin  and  socket  inductances,  which  were  not 
accounted  for  in  the  model.  Likewise,  charge  control  ef- 
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Fig.  7.  (a)  Configuration  used  for  simulation  of  transients  on  a  slripline 
structure  with  nonlinear  loads,  (b)  Model  used  for  AS04  driver  with 
voltage  dependence  of  impedance,  C, -8  pF.  /.(I)  has  the  following 
characteristics:  width  - 108  ns,  rise  time  ••  5  ns,  fall  time  >•  4  ns. 
(c)  Model  used  for  the  input  of  each  AS240  receiver.  The  diode  has  a 
saturation  current  f,  “10”'^  A. 

feels,  which  determine  rise  and  fall  time  degradation  as 
well  as  time  delay  through  the  inverters,  were  not  modeled. 
The  otherwise  overall  good  agreement  indicated  the  valid¬ 
ity  of  the  simplified  device  models. 

Losses  were  also  analyzed  using  (27).  First,  the 
frequency-domain  scattering  parameters  were  calculated 
for  the  transmission  line  shown  in  Fig,  9.  Skin  effect  in  the 
conductor  was  accounted  by  a  ^  behavior  of  the  resis¬ 
tance  per  unit  length.  Once  the  frequency  dependence  of 
the  scattering  parameter  was  determined,  a  fast  Fourier 
transform  (FFT^  algorithm  was  used  to  solve  for  the 
time-domain  Green’s  functions  associated  with  the  S 
parameters.  Then  (20)-(28)  were  used  to  calculate  the 
voltage  waves.  Fig.  9  shows  the  configuration  used  to 
predict  and  compare  waveforms  on  lossless  and  lossy 
microstrip  lines.  The  low-frequency  characteristics  of  the 
line  were;  Lq  »  539  nH/m  and  Q  *  39  pF/m.  The  resis- 
tance  per  unit  length  for  the  lossy  case  was 
=■  //  (in  GHz)  1  kQ/m.  Dielectric  losses  were  neglected 


(a) 


(b) 


Fig.  8.  Comparison  of  theoretical  (plots)  and  experimental  (photo¬ 
graphs)  simulations  for  the  stripline  structure  terminated  with  AS 
inverters  of  Fig.  7  at  (a)  the  near  end  and  (b)  the  far  end.  Vertical 
scales  in  experimental  simulations  have  a  probe  attenuation  factor  of 

to. 
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Pig.  9.  Microstrip  configuration  used  to  simulate  the  effects  of  losses. 
Line  length  -  25  inches  (0.635  meien).  Low-frequency  electrical  char- 
acteristics:  Lg  ••  539  nH/m,  Q  -  39  pF/m.  Loss  characteristics:  -• 
yY  in  GHz  I  kO/m  and  Co  -  0  mhos/meter.  Pulse  characteristics: 
magnitude  ••  4  V,  width  ••  20  ns,  rise  and  fall  limes  ns. 
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Fig.  10.  Comparisons  between  responses  for  lossless  and  lossy  cases  for 
the  microstrip  structure  of  Fig.  9  at  (a)  the  near  end  and  (b)  the  far  end. 


(Go  =  0),  which  is  a  good  representation  of  many  intercon¬ 
nections  in  integrated  circuit  design.  The  above  configura¬ 
tion  was  excited  by  a  4  V  pulse  generator  with  an  internal 
impedance  of  50  S2,  and  a  termination  impedance  of  1  kJ2. 
A  comparison  between  lossy  and  lossless  cases  is  shown  in 
Fig.  10.  As  anticipated,  the  introduction  of  losses  led  to 
rise  and  fall  time  degradation  and  waveform  attenuation. 
Experimental  simulations  were  not  available  due  to  the 
lack  of  reliable  measurement  techniques  for  accurately 
determining  the  frequency  dependence  of  microstrip  loss 
parameters.  Such  information  is  essential  in  developing 
models  to  be  used  in  conjunction  with  the  simulations. 
Other  alternatives  would  include  providing  a  set  of  mea¬ 
sured  scattering  parameters  for  an  arbitrary  (lossy  and 


dispersive)  transmission  line  over  a  wide  frequency  range 
(up  to  18  GHz).  Time-domain  Green’s  functions  can  then 
be  numerically  computed  and  used  to  calculate  the  associ¬ 
ated  response. 


VI.  Conclusion.^ 

This  study  explored  some  important  aspects  of  intercon¬ 
nections  for  digital  and  microwave  applications.  The  prob¬ 
lems  of  losses,  dispersion,  and  load  nonlinearities  were 
analyzed.  A  simple  algorithm  was  derived  for  the  simula¬ 
tion  of  an  arbitrary  time-domain  signal  on  a  structure 
having  all  the  above  properties.  The  algorithm  assumed 
that  the  frequency-dependent  characteristics  of  the  line 
were  available  as  well  as  large-signal  models  for  the 
terminations.  Future  work  includes  the  derivation  of  a 
suitable  loss  and  dispersion  model  and  the  extension  of  the 
:«!goiithm  for  modeling  n-Iine  multiconductor  systems. 
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